setwd("/Users/isaiahmireles/Desktop/Misconceptions")
cohorts <- read.csv("cohorts.csv")
cohorts$exam <- factor(cohorts$exam)
library(tidyverse)

1 Data Structure :

  • notice we have a hierarchical structure to our data.

  • Repeated obs. occur at the student level

    • ie. std. j took responded to item i for both the Midterm and Final
# verify N-ct diagram above
cohorts |> 
  group_by(term, exam) |> 
  summarize(std_ct = n())

  • exams vary in std. ct

  • Notice students after taking the Midterm , drop S10

# Std. dropoff (After Data Clean Filter)
cohorts |>
  group_by(term, exam) |>
  # careful w/ summarize preserving grouping 
  summarize(
    std_ct = n(),
    .groups = "drop"
  ) |>
  pivot_wider(
    names_from = exam,
    values_from = std_ct
  ) |>
  mutate(dropoff = Final - Midterm)
  • notice substantial droppoff from F23
# inspect raw data

2 Per Question Analysis

2.1 Percent Per Question : PPQ

  • the following data is the percent of students that got item i correct

    • not to be confused with pct which defines the percent-score obtained by student j
# pivot_longer
# turn col to rows, rows to col
?pivot_longer
# percent per q dat 
PPQ <- 
  cohorts |>
  pivot_longer(
    # select starts w/ Q, then 1 digit 
    cols = matches("^Q\\d+$"),
    # grab q, 
    names_to = "question",
    values_to = "correct"
  ) |>
  group_by(term, exam, question) |>
  # ct std (matches above)
  summarize(
    n_students = n(),
    PPQ = mean(correct)
  ) 
# Data split PPQ 
PPQ_list <- PPQ |> 
  group_by(term, exam) |> 
  # split into seperate data 
  group_split() |>
  # name the data -- vec. input
  setNames(
    PPQ |> 
      distinct(term, exam) |> 
      arrange(term, exam) |>
      # paste distinct strings for naming together
      unite("name", term, exam) |> 
      # grab vector of name 
      pull(name)
  )
?unite

3 unidentified : T, F Class

cohorts |> filter(unidentified==T) |> count(term, exam)
  • notice these are students whom dropped the course after the midterm

3.1 Do they appear to perform systematically different?

cohorts |> 
  filter(term == "F23",exam=="Midterm") |> 
  group_by(unidentified) |>
  summarize(
    n = n(),
    mean = round(mean(Total.Score), 2),
    sd = round(sd(Total.Score), 2)
  )
  • mean and average deviation from mean appears rather similar
cohorts |> 
  filter(term == "F23",exam=="Midterm") |> 
  summarize(
    n = n(),
    mean = round(mean(Total.Score), 2),
    sd = round(sd(Total.Score), 2)
  )
  • overall it appears to be similar to the overall aswell
library(patchwork)
p1 <- 
  cohorts |> 
  filter(term == "F23") |> 
  mutate(unidentified = factor(unidentified)) |> 
  ggplot(aes(x=Total.Score, color = unidentified)) + 
  geom_histogram(binwidth = 1)

p2 <-
  cohorts |> 
  filter(term == "F23") |> 
  mutate(unidentified = factor(unidentified)) |> 
  ggplot(aes(x=Total.Score, color = unidentified)) + 
  geom_boxplot()


p1 / p2

  • These students appear to performs quite similarly
  • appears to be 2 samples of the same underlying dist.

3.2 Statistical Inference

3.2.1 Hypothesis : Difference of means

\[ H_O : \mu_{\text{identified}}-\mu_{\text{UN-identified}}=0 \\ H_A : \mu_{\text{identified}}-\mu_{\text{UN-identified}}\ne0 \]

# subset dat.
f23_mid <- 
  cohorts |>
  filter(term == "F23", exam == "Midterm")

3.2.2 Test Statistic ( t-dist )

Welch Two Sample t-test

Assumptions

# T-tst comparing 2 groups means
t.test(
  Total.Score ~ unidentified, 
  alternative = "two.sided", 
  data = f23_mid)
## 
##  Welch Two Sample t-test
## 
## data:  Total.Score by unidentified
## t = -0.072867, df = 469.05, p-value = 0.9419
## alternative hypothesis: true difference in means between group FALSE and group TRUE is not equal to 0
## 95 percent confidence interval:
##  -0.9769727  0.9071074
## sample estimates:
## mean in group FALSE  mean in group TRUE 
##            27.30303            27.33796
?t.test
  • there does not appear to be a statistically significant difference
library(effectsize)
cohens_d(
  Total.Score ~ unidentified,
  pooled_sd = FALSE,
  data = f23_mid
)
?cohens_d
  • here we see 0 \(\in\) CI-95%

4 Exam, Midterm Relationship :

cohorts <- 
  cohorts |> 
  filter(unidentified!=T) |> 
  select(-unidentified)
  • We will exclude such students so we may map students across exam (Midterm and Final)
cohorts |> 
  group_by(term, exam) |> 
  summarize(avr = mean(Total.Score))
cohorts_wide <- cohorts |>
  filter(exam %in% c("Midterm", "Final")) |>
  select(term, student_id, exam, pct) |>
  pivot_wider(
    names_from = exam,
    values_from = pct
  )

cohorts_wide |>
  ggplot(aes(x = Midterm, y = Final)) +
  geom_point(alpha = 0.5) +
  geom_smooth(method = "lm", se = TRUE) +
  facet_wrap(~ term) +
  labs(
    x = "Midterm Performance",
    y = "Final Performance"
  ) +
  theme_minimal()
## `geom_smooth()` using formula = 'y ~ x'
## Warning: Removed 14 rows containing non-finite outside the scale range
## (`stat_smooth()`).
## Warning: Removed 14 rows containing missing values or values outside the scale range
## (`geom_point()`).

cohorts |>
  ggplot(aes(x = Total.Score)) +
  geom_histogram(binwidth =  1) +
  facet_grid(term ~ exam)

cohorts |>
  ggplot(aes(x = pct)) +
  geom_density() +
  facet_grid(term ~ exam)

  • appears to be bi-modal behavior

    • Perhaps a mixed beta dist

  • underlying latent variable

4.1 Linear Models

lm_by_term <- cohorts_wide |>
  drop_na(Midterm, Final) |>
  group_by(term) |>
  group_nest() |>
  mutate(
    results = map(data, ~ {
      
      fit <- lm(Final ~ Midterm, data = .x)
      
      list(
        model = fit,
        summary = summary(fit)
      )
    })
  ) |>
  select(term, results) |>
  deframe()
lm_by_term$F23$model
## 
## Call:
## lm(formula = Final ~ Midterm, data = .x)
## 
## Coefficients:
## (Intercept)      Midterm  
##      0.1848       0.8196
lm_by_term$W24$model
## 
## Call:
## lm(formula = Final ~ Midterm, data = .x)
## 
## Coefficients:
## (Intercept)      Midterm  
##      0.1323       0.8359
lm_by_term$S24$model
## 
## Call:
## lm(formula = Final ~ Midterm, data = .x)
## 
## Coefficients:
## (Intercept)      Midterm  
##     0.01273      0.95441

5 Mixed Eff. Model

cohorts |>
  group_by(term, student_id) |>
  arrange(term, student_id, factor(exam, levels = c("Midterm", "Final")))
library(lme4)
## Warning: package 'lme4' was built under R version 4.4.3
## Loading required package: Matrix
## 
## Attaching package: 'Matrix'
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
mixed_model <- cohorts |>
  mutate(
    exam = factor(exam, levels = c("Midterm", "Final")),
    term = factor(term)
  ) |>
  lmer(
    pct ~ exam * term + (1 | term:student_id),
    data = _
  )
mixed_model
## Linear mixed model fit by REML ['lmerMod']
## Formula: pct ~ exam * term + (1 | term:student_id)
##    Data: mutate(cohorts, exam = factor(exam, levels = c("Midterm", "Final")),  
##     term = factor(term))
## REML criterion at convergence: -1723.063
## Random effects:
##  Groups          Name        Std.Dev.
##  term:student_id (Intercept) 0.14429 
##  Residual                    0.08039 
## Number of obs: 1474, groups:  term:student_id, 744
## Fixed Effects:
##       (Intercept)          examFinal            termS24            termW24  
##           0.80285            0.04009           -0.02848           -0.03671  
## examFinal:termS24  examFinal:termW24  
##          -0.06213           -0.03325
summary(mixed_model)
## Linear mixed model fit by REML ['lmerMod']
## Formula: pct ~ exam * term + (1 | term:student_id)
##    Data: mutate(cohorts, exam = factor(exam, levels = c("Midterm", "Final")),  
##     term = factor(term))
## 
## REML criterion at convergence: -1723.1
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.99108 -0.39591  0.09085  0.45928  2.32793 
## 
## Random effects:
##  Groups          Name        Variance Std.Dev.
##  term:student_id (Intercept) 0.020819 0.14429 
##  Residual                    0.006463 0.08039 
## Number of obs: 1474, groups:  term:student_id, 744
## 
## Fixed effects:
##                    Estimate Std. Error t value
## (Intercept)        0.802854   0.009575  83.850
## examFinal          0.040092   0.006635   6.042
## termS24           -0.028475   0.016576  -1.718
## termW24           -0.036707   0.013548  -2.710
## examFinal:termS24 -0.062125   0.011504  -5.400
## examFinal:termW24 -0.033250   0.009399  -3.538
## 
## Correlation of Fixed Effects:
##             (Intr) exmFnl trmS24 trmW24 eF:S24
## examFinal   -0.343                            
## termS24     -0.578  0.198                     
## termW24     -0.707  0.243  0.408              
## exmFnl:tS24  0.198 -0.577 -0.342 -0.140       
## exmFnl:tW24  0.242 -0.706 -0.140 -0.342  0.407
library(performance)
## Warning: package 'performance' was built under R version 4.4.3
icc(mixed_model)

Null Model :

null_model <- cohorts |>
  mutate(
    term = factor(term)
  ) |>
  lmer(
    pct ~ 1 + (1 | term:student_id),
    data = _
  )

performance::icc(null_model)

6 Exploratory Factor Analysis

6.1 Scree Plot*

library(tidyverse)
library(psych)
## 
## Attaching package: 'psych'
## The following object is masked from 'package:effectsize':
## 
##     phi
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
scree_data <- cohorts |>
  group_by(term, exam) |>
  group_modify(~ {

    items <- .x |>
      select(starts_with("Q")) |>
      # Remove items that were not administered for this exam
      select(where(~ !all(is.na(.x)))) |>
      # Remove items with no variation
      select(where(~ n_distinct(.x, na.rm = TRUE) > 1))

    R <- psych::tetrachoric(items)$rho

    eig <- eigen(R)$values

    tibble(
      component = seq_along(eig),
      eigenvalue = eig
    )
  }) |>
  ungroup()
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
## Warning in cor.smooth(mat): Matrix was not positive definite, smoothing was
## done
scree_data |>
  ggplot(aes(x = component, y = eigenvalue)) +
  geom_point() +
  geom_line() +
  geom_hline(yintercept = 1, linetype = "dashed") +
  facet_grid(term ~ exam, scales = "free_x") +
  labs(
    x = "Component",
    y = "Eigenvalue",
    title = "Scree plots by term and exam"
  )

  • for each exam

6.2 IRT Modeling