Machine Learning/HMM R Example

Examples of using HMM R packages, based on the model in "A Bayes Net Toolkit for Student Modeling in Intelligent Tutoring Systems" by Chang, et. al. We're trying to come up with an estimate for how well each student knows a certain area of knowledge (which we're calling a skill). We observe each student's performance on answering some number of questions that use this skill, and mark whether they got them correct or incorrect.

We assume that at each time point, a student is in one of two states: either they "know" the skill, or they "do not know" the skill. If they know they skill, they are more likely to generate a correct output; if not, they are less likely; but in each case, it is stochastic (a student has a probability of guessing the correct answer even if they don't know the skill, and of slipping/getting it wrong even if they do know the skill). Between each time point, there is a transition probability from know -> don't know (forgetting, which Change et. al constrain to 0) and from don't know -> know (learning). Finally, there is a probability that the student enters already knowing the skill. So we have five parameters: two transition probabilities (learn and forget), two outcome probabilities based on state (guess and slip), and initial state probabilities (already know).

The data (student_outcomes.csv) is for a single skill, measuring various students' performance on that skill: a series of correct/incorrect responses, at various times. We're ignoring the time data for the moment (other than for ordering purposes), and trying to fit the HMM model. Once we have it, we can then figure out, for each student, an estimated likelihood of being in the "know" state at their last observed output.

hmm.discnp

require("hmm.discnp")
student_outcomes = read.csv("student_outcomes.csv", header=TRUE)

# convert created_at from a string
student_outcomes$created_at = as.POSIXct(as.character(student_outcomes$created_at))

# remove users with few observations on this skill
by_user = split(student_outcomes, student_outcomes$student_id)
obs_by_user = sapply(by_user, nrow)
valid_users = names(obs_by_user[obs_by_user > 10])
student_outcomes = student_outcomes[student_outcomes$student_id %in% valid_users,]

by_good_user = split(student_outcomes, student_outcomes$student_id)

# attempt to estimate model parameters
my_hmm = hmm(by_good_user, yval=c(0,1),
    par0=list(tpm=rbind(c(0.8,0.2),c(0.01,0.99)),
              Rho=rbind(c(0.75,0.25),c(0.25,0.75))),
    stationary=FALSE)
if (!my_hmm$converged) {
  print(sprintf("Error!  HMM did not converge for skill %s!", skill))
} else {
  for (user_id in valid_users) {
    student_est = sp(correct_by_user[[user_id]], object = my_hmm, means=TRUE)
    print(sprintf("%s/%s: %f chance know, %f chance correct", skill, user_id, student_est$probs[2,ncol(student_est$probs)], student_est$means[length(student_est$means)]))
    # print(correct_by_user[[user_id]])
  }
}
# transition probability matrix
my_hmm$tpm
# output probabilities
my_hmm$Rho
# initial probabilities (don't know/know)
my_hmm$ispd

msm

student_outcomes = read.csv("student_outcomes.csv", header=TRUE)

# convert created_at from a string
student_outcomes$created_at = as.POSIXct(as.character(student_outcomes$created_at))

# remove users with few observations on this skill
min_observations = 10
by_user = split(student_outcomes, student_outcomes$student_id)
obs_by_user = sapply(by_user, nrow)
valid_users = names(obs_by_user[obs_by_user >= min_observations])
student_outcomes = student_outcomes[student_outcomes$student_id %in% valid_users,]

require("msm")
# convert time to simple sequence
student_outcomes$created_index = c(sapply(by_user, function(df) {1:nrow(df)}), recursive=TRUE)
my_hmm = msm(correct ~ created_index, subject = student_id, data = student_outcomes,
                qmatrix = rbind(c(NA,0.25),c(0.25,NA)),
                hmodel = list(hmmBinom(1,0.3), hmmBinom(1,0.7)),
                obstype = 2,
                initprobs = c(0.5,0.5),
                est.initprobs = TRUE,
                method="BFGS"
               )
# display final probability for each user
for (user_id in valid_users) {
    student_est = estimate_knowledge(correct_by_user[[user_id]], my_msm)
    print(sprintf("%s/%s: %f chance know, %f chance correct", skill, user_id, student_est[["p_know"]], student_est[["p_correct"]]))
    print(correct_by_user[[user_id]])
}