1 Preliminary

library(tidyverse)
rm(list=ls()); setwd("/Users/isaiahmireles/Desktop/Misconceptions"); load("S10 Data.RData")
cohorts <- cohorts |> rename(rw_num = X) 
# de-identify year due to such recency  
cohorts <- 
  cohorts |>
  mutate(
    recency = case_when(
      stringr::str_detect(as.character(term), "S") ~ 1L,
      stringr::str_detect(as.character(term), "W") ~ 2L,
      stringr::str_detect(as.character(term), "F") ~ 3L,
      .default = NA_integer_
    ),
    term = case_when(
      stringr::str_detect(as.character(term), "F") ~ "Fall",
      stringr::str_detect(as.character(term), "W") ~ "Winter",
      stringr::str_detect(as.character(term), "S") ~ "Spring",
      .default = NA_character_
    )
  )
summary_tbl <- 
  cohorts |>
  filter(exam == "Midterm") |>
  group_by(term) |>
  summarize(
    xbar = mean(Total.Score),
    xbar_pct = mean(pct_rounded),
    med = median(Total.Score),
    med_pct = median(pct_rounded),
    sd_score = sd(Total.Score),
    sd_pct = sd(pct_rounded)
  ) |>
  mutate(across(where(is.numeric), ~ round(.x, 3)))

summary_tbl
Midterm_cohorts <- 
  cohorts |>
  filter(exam == "Midterm") |>
  group_by(term) |>
  mutate(
    below_median = pct_rounded < median(pct_rounded)
  ) |>
  ungroup()

Now, I add a feature based on the following : D+ or below (.65)

# fix this feature later
Midterm_cohorts <- Midterm_cohorts |> mutate(D_or_below = pct_rounded <= .65)
  • purpose : study class characteristic(s)

    • Ie. are the behavior(s) [Correct/Incorrect-ness patterns] different?

2 Analysis

2.1 Tetrachoric correlation(s)

  • here we have binary data, hence NOT using pearson-correlation

Ai Generated Code* :

library(psych)
# Find column names like Q1, Q2, ..., Q34
question_cols <- grep(
  "^Q\\d+$",
  names(Midterm_cohorts),
  value = TRUE
)

# Split the Midterm records by term and calculate one matrix per term
tetra_by_term <- split(
  Midterm_cohorts,
  Midterm_cohorts$term
) |>
  lapply(function(dat) {
    # Keep only the item responses
    items <- dat |>
      select(all_of(question_cols)) |>
      as.data.frame()

    # Keep items with at least two observed response values
    has_variation <- vapply(
      items,
      function(x) length(unique(stats::na.omit(x))) > 1,
      logical(1)
    )

    # Remove all-missing or constant items
    items <- items[, has_variation, drop = FALSE]

    # Estimate tetrachoric correlations among the binary items
    tetrachoric(
      items,
      correct = 0.5, # Correct for item pairs with a zero-count cell
      smooth = TRUE, # Adjust the matrix if it is not positive definite
      global = TRUE  # Use shared threshold estimates across item pairs
    )$rho
  })
## 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
  • “positive definite” issue
# Check the number of items included in each term's matrix
lapply(tetra_by_term, dim) # nice 
## $Fall
## [1] 34 34
## 
## $Spring
## [1] 34 34
## 
## $Winter
## [1] 33 33

2.1.1 Saving PCs

Ai Generated Code*

# Run principal-components analysis on each term's tetrachoric matrix
# lapply() preserves the names from tetra_by_term in the resulting list
pc <- lapply(tetra_by_term, function(R) {
  stats::princomp(covmat = R)
})

# Optional: save the full named list to reload in a later R session
# saveRDS(pc, "pc_by_term.rds")
# ct PCs > 1
pcs_above_1 <- data.frame(
  term = names(pc),
  n_pcs_above_1 = vapply(
    pc,
    function(pca) sum(pca$sdev^2 > 1),
    numeric(1)
  ),
  row.names = NULL
)

# View the count for each term
pcs_above_1

2.2 Scree Plot(s)

# Load ggplot2 to create reusable plot objects
library(ggplot2)

# Build a scree plot object from one term's tetrachoric matrix
make_scree_plot <- function(term) {
  # Convert the correlation matrix to a principal-components object
  pc <- stats::princomp(covmat = tetra_by_term[[term]])

  # Square the component standard deviations to get eigenvalues
  plot_data <- data.frame(
    component = seq_along(pc$sdev),
    eigenvalue = pc$sdev^2
  )

  # Return a ggplot object for this term
  ggplot(plot_data, aes(x = component, y = eigenvalue)) +
    geom_line() +
    geom_point() +
    geom_hline(yintercept = 1, linetype = "dashed", color = "red") +
    scale_x_continuous(breaks = plot_data$component) +
    labs(
      title = paste("Scree Plot —", term, "Midterm"),
      x = "Component number",
      y = "Eigenvalue"
    ) +
    theme_minimal()
}

# Create one named plot object for each term
Scree_Fall <- make_scree_plot("Fall")
Scree_Winter <- make_scree_plot("Winter")
Scree_Spring <- make_scree_plot("Spring")
# Save each plot object so it can be loaded in a later R session
saveRDS(Scree_Fall, "Scree_Fall.rds")
saveRDS(Scree_Winter, "Scree_Winter.rds")
saveRDS(Scree_Spring, "Scree_Spring.rds")
  • alternatively, we may use screeplot()

2.2.1 Global

library(patchwork)
Scree_Fall / Scree_Winter / Scree_Spring

  • globally, they appear to hold a similar trend across all exam with spring being a little particular

2.2.2 Fall

Scree_Fall

  • Qty above \(\lambda_i > 1\) (Fall): 9
    • At most many factors^
  • Based on elbow, it looks like maybe 2-3 factors
  • One dominant dimension

2.2.3 Winter

Scree_Winter

  • Qty above \(\lambda_i > 1\) (Winter): 11
    • At most many factors^
  • Based on elbow, it looks like maybe 2-3 factors
  • One dominant dimension

2.2.4 Spring

Scree_Spring

  • Qty above \(\lambda_i > 1\) (Spring): 12
    • At most many factors^
  • Based on elbow, it looks like maybe 2-3 factors
  • One dominant dimension

2.2.4.1 Summary

  • it looks like realistically, maybe there are 2 factors (at most); based on elbow; but realistically it looks like based on elbow, just one dominant dimension.

2.3 Factor Loading on particular Problem(s) :

Ai Generated Code*

# Extract every item loading from each term's PCA and combine them
pc_loadings_by_term <- bind_rows(lapply(names(pc), function(term) {
  # Convert this term's loadings to a regular matrix
  loading_matrix <- as.matrix(pc[[term]]$loadings)

  # Reshape the matrix so each row gives one item's loading on one component
  as.data.frame(as.table(loading_matrix)) |>
    rename(
      item = Var1,
      component = Var2,
      loading = Freq
    ) |>
    mutate(
      term = term,
      abs_loading = abs(loading)
    )
}))

# For each term and component, keep the five items with the largest
# absolute loadings; keep the signed loading for interpretation
top_items_by_component <- pc_loadings_by_term |>
  group_by(term, component) |>
  slice_max(
    order_by = abs_loading,
    n = 5,
    with_ties = FALSE
  ) |>
  arrange(term, component, desc(abs_loading))
# View the strongest-loading items for each component in each term
top_items <- 
  top_items_by_component |>
  filter(component %in% c("Comp.1", "Comp.2")) |>
  ungroup() |>
  arrange(term, component, desc(abs_loading))

top_items

Aside :

top_items_by_component |> 
  filter(component=="Comp.1"|component=="Comp.2") |> 
  group_by(term) |> 
  arrange(abs_loading) |> 
  count(term)
  • Here we see 10 obs. per term for either Comp.1 or, Comp.2
top_items_by_component |> 
  filter(component=="Comp.1"|component=="Comp.2") |> 
  group_by(term) |> 
  arrange(abs_loading) |> 
  count(term, component)
  • whereby we see 5 for (term, component) pairs

    • ie. each comp has about 5 obs

2.3.1 Are these similar questions?

Recall :

top_items |> select(term, item, component)
  • here I just care about which factors & questions, not the strength of their loadings.
# Combine the three tables into one and standardize the join columns
midterm_details <- bind_rows(Midterms_lst) |>
  mutate(
    term = case_when(
      as.character(term) == "F23" ~ "Fall",
      as.character(term) == "W24" ~ "Winter",
      as.character(term) == "S24" ~ "Spring",
      .default = as.character(term)
    ),
    question = as.character(question)
  )
# Add Midterm details to the three selected columns
top_items_details <- top_items |>
  ungroup() |>
  select(term, item, component) |>
  mutate(item = as.character(item)) |>
  left_join(
    midterm_details,
    by = c("term" = "term", "item" = "question"),
    relationship = "many-to-one"
  )

# Show every row and column
top_items_details |> 
  select(
    term, item, component, 
    Question_Order, PPQ, difference, 
    description, comment)
  • here i have joined my Midterms_lst details to be able to compare

    • Apples to apples and not oranges to oranges

We will take a closer look during actual meeting due to the confidentiality of exam materials.

2.4 Parallel Analysis :

  • i dont know what this is

2.5 IRT Modeling