| Title: | Random Effects Latent Class Analysis |
|---|---|
| Description: | Fits standard and random effects latent class models. The single level random effects model is described in Qu et al <doi:10.2307/2533043> and the two level random effects model in Beath and Heller <doi:10.1177/1471082X0800900302>. Examples are given for their use in diagnostic testing. |
| Authors: | Ken Beath [aut, cre] |
| Maintainer: | Ken Beath <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 1.1-4 |
| Built: | 2024-12-23 06:46:47 UTC |
| Source: | CRAN |
Returns AIC for a randomLCA object.
## S3 method for class 'randomLCA' AIC(object, ..., k = 2)## S3 method for class 'randomLCA' AIC(object, ..., k = 2)
object |
randomLCA object |
... |
additional argument; currently none is used. |
k |
penalty per parameter |
AIC
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(AIC(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(AIC(pap.lca2))
Returns AIC with penalty 3 for a randomLCA object.
AIC3(object)AIC3(object)
object |
randomLCA object |
AIC3.
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(AIC3(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(AIC3(pap.lca2))
Returns BIC for a randomLCA object.
## S3 method for class 'randomLCA' BIC(object, ...)## S3 method for class 'randomLCA' BIC(object, ...)
object |
randomLCA object |
... |
additional argument; currently none is used. |
BIC
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(BIC(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(BIC(pap.lca2))
The conditional probabilities are obtained integrating over the period random effect.
calcCond2Prob(object, conditionalp = 0.5)calcCond2Prob(object, conditionalp = 0.5)
object |
RandomLCA object |
conditionalp |
the percentiles for the random effect |
Returns a data frame containing class, block, outcome, outcomep (outcome probability) and perc (percentiles of the random effect) if conditionalp is specified. For example a conditionalp of 0.5 is the 50th percentile or the median corresponding to a random effect of zero. 0.025 and 0.975 correspond to the 2.5th and 97.5th percential, so the region between them if 95% of the variation in the data.
Ken Beath [email protected]
symptoms.lca2random2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, level2 = TRUE, nclass = 2, level2size = 4, constload = FALSE, cores = 1) print(calcCond2Prob(symptoms.lca2random2))symptoms.lca2random2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, level2 = TRUE, nclass = 2, level2size = 4, constload = FALSE, cores = 1) print(calcCond2Prob(symptoms.lca2random2))
Calculates the conditional outcome probabilities for random effects models or for standard latent class returns the outcome probabilities. For random effects, the outcome probabilities may be calculated for various percentiles of the random effect.
calcCondProb(object, conditionalp = 0.5)calcCondProb(object, conditionalp = 0.5)
object |
RandomLCA object |
conditionalp |
the percentiles for the random effect |
Returns a data frame containing class, block, outcome, outcomep (outcome probability) and perc (percentiles of the random effect) if conditionalp is specified. For example a conditionalp of 0.5 is the 50th percentile or the median corresponding to a random effect of zero. 0.025 and 0.975 correspond to the 2.5th and 97.5th percential, so the region between them if 95% of the variation in the data.
Ken Beath [email protected]
dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(calcCondProb(dentistry.lcarandom))dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(calcCondProb(dentistry.lcarandom))
Calculates the marginal outcome probabilities for a random effects latent class model, by integrating the outcome probability over the random effect. This is performed using Gauss-Hermite quadrature with the number of quadrature points specified for the model fitting.
calcMargProb(object)calcMargProb(object)
object |
randomLCA object |
Returns a data frame containing class, block, outcome, outcomep (outcome probability).
Ken Beath
dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(calcMargProb(dentistry.lcarandom))dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(calcMargProb(dentistry.lcarandom))
The class probabililities for the model are returned.
classProbs(object)classProbs(object)
object |
randomLCA object |
Simply extracts the corresponding variable from the randomLCA object.
A vector of class probabilities for each class.
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(classProbs(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(classProbs(pap.lca2))
Six dentists evaluated dental x-rays for incipient caries in Handelman et al (1986), data consistss of 5 of the dentists analysed by Espeland and Handelman (1989) using a latent class model. Further analysis incorprating a random effects latent class model by Qu et al (1996), and by Albert and Dodd (2004)
dentistrydentistry
A data frame with 32 observations on the following 6 variables.
V1Dentist 1
V2Dentist 2
V3Dentist 3
V4Dentist 4
V5Dentist 5
freqNumber of subjects
Espeland and Handelman (1989)
Handelman, S.L., Leverett, D.H., Espeland, M.A. and Curzon, J.A. (1986) Clinical radiographic evaluation of sealed carious and sound tooth surfaces. Journal of the American Dental Association, 113, 751–754.
Espeland, M.A. and Handelman, S.L. (1989) Using latent class models to characterize and assess relative error in discrete distributions. Biometrics, 45, 587–599.
Qu, Y., Tan, M. and Kutner, M.H. (1996) Random effects models in latent class analysis for evaluating accuracy of diagnostic tests. Biometrics, 52, 797–810.
Albert P.S. and Dodd, L.E. (2004) A cautionary note on the robustness of Latent Class Models for estimating diagnostic error without a gold standard. Biometrics, 60, 427–435.
# fit LCR model from Qu et al (1996) dentistry.lca <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, cores = 1) # start with constant laoding dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) # allow loading to vary by dentist dentistry.lcarandomunequal <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, constload = FALSE, probit = TRUE, cores = 1)# fit LCR model from Qu et al (1996) dentistry.lca <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, cores = 1) # start with constant laoding dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) # allow loading to vary by dentist dentistry.lcarandomunequal <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, constload = FALSE, probit = TRUE, cores = 1)
Extract fitted values for randomLCA object.
## S3 method for class 'randomLCA' fitted(object, ...)## S3 method for class 'randomLCA' fitted(object, ...)
object |
randomLCA object |
... |
additional argument; currently none is used. |
A data frame. The first columns of the data frame correspond to the patterns, followed by the frequency of each pattern, and then the fitted number for each pattern.
Ken Beath [email protected]
dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(fitted(dentistry.lcarandom))dentistry.lcarandom <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, random = TRUE, probit = TRUE, cores = 1) print(fitted(dentistry.lcarandom))
Opinions collected on gender roles in a study by Felling et al (1987). This was originally published in Heinen (1996) and subsequently in Galindo Garre and Vermunt (2006).
genderrolegenderrole
A data frame with 16 observations on the following 5 variables.
Q1Women's liberation sets women against men.
Q2It's better for a wife not to have a job because that always poses problems in the household, especially if there are children.
Q3The most natural situation occurs when the man is the breadwinner and the woman runs the household and takes care of the children.
Q4It isn't really as important for a girl to get a good education as it is for a boy.
Q5A woman is better suited to raise small children than a man.
FreqNumber of subjects
Galindo Garre and Vermunt (2006)
Felling, A., Peters, J., and Schreuder, O. (1987) Religion in Dutch society 85: Documentation of a national survey on religious and secular attitudes in 1985. Amsterdam: Steinmetz Archive.
Galindo Garre, F. and Vermunt, J.K. (2006) Avoiding boundary estimates in latent class analysis by Bayesian posterior mode estimation. Behaviormetrika, 33, 43–59.
Heinen, T. (1996) Latent Class and Discrete Latent Trait Models: Similarities and Differences.
# standard latent class genderrole.lca1 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, cores = 1) genderrole.lca2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, cores = 1) genderrole.lca3 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, cores = 1) # repeat with random effect with constant loading # increase quadrature points and/or use higher penalty to obtain # convergence genderrole.lca1random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, random = TRUE, cores = 1) genderrole.lca2random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, random = TRUE, penalty = 0.1, quadpoints = 61, cores = 1) genderrole.lca3random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, random = TRUE, penalty = 0.1, quadpoints = 61, cores = 1) # improved BIC for 1 class random print(c(BIC(genderrole.lca1), BIC(genderrole.lca2), BIC(genderrole.lca3))) print(c(BIC(genderrole.lca1random), BIC(genderrole.lca2random), BIC(genderrole.lca3random))) # can also repeat fits without constant loading to give mixture of IRT models genderrole.lca1random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, random = TRUE, constload = FALSE, cores = 1) genderrole.lca2random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, random = TRUE, constload = FALSE, quadpoints = 61, penalty = 0.1, cores = 1) genderrole.lca3random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, random = TRUE, constload = FALSE, quadpoints = 61, penalty = 0.1, cores = 1) # no improvement in fit print(c(BIC(genderrole.lca1random2), BIC(genderrole.lca2random2), BIC(genderrole.lca3random2)))# standard latent class genderrole.lca1 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, cores = 1) genderrole.lca2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, cores = 1) genderrole.lca3 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, cores = 1) # repeat with random effect with constant loading # increase quadrature points and/or use higher penalty to obtain # convergence genderrole.lca1random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, random = TRUE, cores = 1) genderrole.lca2random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, random = TRUE, penalty = 0.1, quadpoints = 61, cores = 1) genderrole.lca3random <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, random = TRUE, penalty = 0.1, quadpoints = 61, cores = 1) # improved BIC for 1 class random print(c(BIC(genderrole.lca1), BIC(genderrole.lca2), BIC(genderrole.lca3))) print(c(BIC(genderrole.lca1random), BIC(genderrole.lca2random), BIC(genderrole.lca3random))) # can also repeat fits without constant loading to give mixture of IRT models genderrole.lca1random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 1, random = TRUE, constload = FALSE, cores = 1) genderrole.lca2random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, random = TRUE, constload = FALSE, quadpoints = 61, penalty = 0.1, cores = 1) genderrole.lca3random2 <- randomLCA(genderrole[, 1:5], freq = genderrole$Freq, nclass = 3, random = TRUE, constload = FALSE, quadpoints = 61, penalty = 0.1, cores = 1) # no improvement in fit print(c(BIC(genderrole.lca1random2), BIC(genderrole.lca2random2), BIC(genderrole.lca3random2)))
Serum samples are tested for HIV by 4 different biossays in Alvord et al (1988) and sensitivity and specificity determined using latent class analysis. Qu et al (1996) repeat the analysis using a model incorporating a random effect.
hivtestshivtests
A data frame with 16 observations on the following 5 variables.
V1Test 1
V2Test 2
V3Test 3
V4Test 4
freqNumber of subjects
Qu, Tan and Kutner (1989)
Alvord, W.G., Drummond, J.E., Arthur, L.O., Goedert, J.J., Levine, P.H., Murphy, E.L., Weiss, S.H., and Blattner, W.A. (1988) A method for predicting individual HIV infection status in the absence of clinical information. AIDS Research and Human Retroviruses, 4, 295–304.
Qu, Y., Tan, M. and Kutner, M.H. (1996) Random effects models in latent class analysis for evaluating accuracy of diagnostic tests. Biometrics, 52, 797–810.
# fit standard latent class hivtests.lca2 <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, cores = 1) # with random effect and constant loading hivtests.lca2random <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, random = TRUE, quadpoints = 101, penalty = 1.0, cores = 1) # with random effect and variable loading # for this model there are 13 parameters fitted to 16 observations, so model is fairly unstable hivtests.lca2random2 <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, random = TRUE, constload = FALSE, quadpoints = 101, penalty = 1.0, cores = 1) # BIC shows best model is random effects with constant loading print(c(BIC(hivtests.lca2), BIC(hivtests.lca2random), BIC(hivtests.lca2random2)))# fit standard latent class hivtests.lca2 <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, cores = 1) # with random effect and constant loading hivtests.lca2random <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, random = TRUE, quadpoints = 101, penalty = 1.0, cores = 1) # with random effect and variable loading # for this model there are 13 parameters fitted to 16 observations, so model is fairly unstable hivtests.lca2random2 <- randomLCA(hivtests[, 1:4], freq = hivtests$freq, random = TRUE, constload = FALSE, quadpoints = 101, penalty = 1.0, cores = 1) # BIC shows best model is random effects with constant loading print(c(BIC(hivtests.lca2), BIC(hivtests.lca2random), BIC(hivtests.lca2random2)))
Returns log Likelihood for a randomLCA object.
## S3 method for class 'randomLCA' logLik(object, ...)## S3 method for class 'randomLCA' logLik(object, ...)
object |
randomLCA object |
... |
additional argument; currently none is used. |
The loglikelihood.
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(logLik(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(logLik(pap.lca2))
For each observation the posterior class probability is determined for each class, and then the class with the maximum posterior class probability is returned.
maxPostClass(object)maxPostClass(object)
object |
randomLCA object |
Returns the class with the maximum posterior class probability for each observation.
A data frame. The first columns of the data frame correspond to the patterns, followed by the frequency of each pattern, and then the class with the maximum posterior class probability. The returned result is for the summarised data. If raw data is used, that is no frequencies, and it is required to calculated the posterior class probability for each observation then it is simply required to merge the maximum class with the raw data, possibly removing any variable "Freq" in the raw data.
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) maxClass.lca2 <- maxPostClass(pap.lca2) names(maxClass.lca2)[length(names(maxClass.lca2))] <- "maxProb.lca2" pap.lca3 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, cores = 1) maxClass.lca3 <- maxPostClass(pap.lca3) names(maxClass.lca3)[length(names(maxClass.lca3))] <- "maxProb.lca3" maxClass <- merge(maxClass.lca2, maxClass.lca3) # aggregate because there is more than one record for each combination # shows the relationship between the 2 and 3 class models print(aggregate(maxClass$Freq, list(maxProb2 = maxClass$maxProb.lca2, maxProb3 = maxClass$maxProb.lca3), sum))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) maxClass.lca2 <- maxPostClass(pap.lca2) names(maxClass.lca2)[length(names(maxClass.lca2))] <- "maxProb.lca2" pap.lca3 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, cores = 1) maxClass.lca3 <- maxPostClass(pap.lca3) names(maxClass.lca3)[length(names(maxClass.lca3))] <- "maxProb.lca3" maxClass <- merge(maxClass.lca2, maxClass.lca3) # aggregate because there is more than one record for each combination # shows the relationship between the 2 and 3 class models print(aggregate(maxClass$Freq, list(maxProb2 = maxClass$maxProb.lca2, maxProb3 = maxClass$maxProb.lca3), sum))
Four tests were performed on hospital patients to determine if a myocardial infarction had occurred.
myocardialmyocardial
A data frame with 32 observations on the following 6 variables.
Q.waveresult from ECG test
Historyclinical history
LDHflipped, enzyme related to tissue breakdown
CPKhigh, creatine kinase or creatine phosphokinase, related to muscle damage
freqNumber of subjects
Rindskopf and Rindskopf (1986)
Galen, R.S. and Gambino, S.R. (1975) Beyond Normality: The Predictive Value and Efficiency of Medical Diagnosis. Wiley:New York.
Rindskopf, D. and Rindskopf, W. (1986) The Value of Latent Class Analysis in Medical Diagnosis. Statistics in Medicine, 5, 21–27.
# fit 2 class model from Rindskopf and Rindskopf (1986) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1)# fit 2 class model from Rindskopf and Rindskopf (1986) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1)
Extract outcome probabilities and confidence intervals for a randomLCA object.
## S3 method for class 'randomLCA' outcomeProbs(object, level = 0.95, boot=FALSE, type="perc", R=999, scale=c("prob","raw"), cores = max(detectCores() %/% 2, 1), ...)## S3 method for class 'randomLCA' outcomeProbs(object, level = 0.95, boot=FALSE, type="perc", R=999, scale=c("prob","raw"), cores = max(detectCores() %/% 2, 1), ...)
object |
randomLCA object |
level |
confidence interval |
boot |
use parametric bootstrap to obtain confidence interval |
type |
type of bootstrap confidence intervals to use, with "perc" or "norm" valid, see boot.ci for description. |
R |
replications for parametric bootstrap |
scale |
either "prob" where probabilities are returned, the default, or "raw" where the probabilities are returned on the logit or probit scale, depending on which scale was selected in the randomLCA function |
cores |
number of cores to use when bootstrapping, shuld be at least 1 less than available cores |
... |
additional argument; currently none is used. |
Confidence intervals are calculated based on asymptotic normality of the estimates transformed by either the inverse of the probit or logistic, or using parametric bootstrap. The asymptotic confidence intervals are currently only available for models without random effects. For the confidence intervals obtained from the parametric bootstrap, the bootstrap is performed on the data that has been transformed to the logit or probit scale, as appropriate.
Data frame consisting of outcome probabilities and confidence intervals. One for each class.
Ken Beath
# standard latent class with 2 classes dentistry.lca2 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, cores = 1) print(outcomeProbs(dentistry.lca2)) # print on the default logit scale print(outcomeProbs(dentistry.lca2, scale = "raw")) # convert back to probabilities print(1.0/(1.0+exp(-outcomeProbs(dentistry.lca2, scale = "raw")[[1]]))) print(1.0/(1.0+exp(-outcomeProbs(dentistry.lca2, scale = "raw")[[2]])))# standard latent class with 2 classes dentistry.lca2 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, cores = 1) print(outcomeProbs(dentistry.lca2)) # print on the default logit scale print(outcomeProbs(dentistry.lca2, scale = "raw")) # convert back to probabilities print(1.0/(1.0+exp(-outcomeProbs(dentistry.lca2, scale = "raw")[[1]]))) print(1.0/(1.0+exp(-outcomeProbs(dentistry.lca2, scale = "raw")[[2]])))
The Positive Action program is a series of interventions designed to reduce negative behaviours in elementary-school students. In a study in Hawaii (Beets et al, 2006) information was recorded from students in the treatment group about whether the various parts of the program were implemented. While it is useful to describe the proportion of students experiencing implementation of each part of the program, latent class analysis will reveal if there are specific patterns to the implementation of the program (Alcock, 2008).
pappap
A data frame with 606 observations summarising the answers for 1566 students on the following 11 variables. For each variable the data has been dichotomized so that a 0 represents no implementation and a 1 represents some implementation.
Q1you receive stickers from your teacher for doing positive actions?
Q2you receive a word of the week card from your teacher?
Q3you put notes in an icu box?
Q4your teacher read notes about you from the icu box?
Q5your teacher read your notes from the icu box?
Q6your class receive a token for meeting your classroom goals?
Q7you participate in a positive action assembly?
Q8your class receive a balloon in an assembly for achieving their classroom goals?
Q9your class participate in whole school positive action celebrations?
Q10most weeks were you taught a positive action lesson?
FreqNumber of subjects
Alcock (2008)
Alcock, A. (2008). Latent Class Analysis (LCA) and Latent Profile Analysis (LPA). Retrieved from http://people.oregonstate.edu/~acock/growth/Mplus files/Latent Class Analysis Presentation.doc
Beets, M. W., Flay, B. R., Vuchinich, S., Snyder, F. J., Acock, A., Li, K. K., Durlak, J. (2009). Use of a social and character development program to prevent substance use, violent behaviors, and sexual activity among elementary-school students in Hawaii. American Journal of Public Health, 99(8), 1438-1445
# standard latent class pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) # standard latent class pap.lca3 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, cores = 1) pap.lca4 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 4, cores = 1) pap.lca5 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 5, cores = 1) # repeat with random effect with constant loading # once BIC increases fitting further models is unnecessary pap.lca1random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, random = TRUE, cores = 1) pap.lca2random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, random = TRUE, cores = 1) pap.lca3random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, random = TRUE, cores = 1) # can also repeat fits without constant loading to give mixture of IRT models pap.lca1random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, random = TRUE, constload = FALSE, cores = 1) pap.lca2random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, random = TRUE, constload = FALSE, cores = 1) pap.lca3random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, random = TRUE, constload = FALSE, cores = 1) # produce table of BIC values # shows 4 class best of standard latent class # but 2 class latent class with constant loading has better BIC pap.bic <- data.frame(bic = c(BIC(pap.lca1), BIC(pap.lca2), BIC(pap.lca3), BIC(pap.lca4), BIC(pap.lca5)), bic2 = c(BIC(pap.lca1random), BIC(pap.lca2random), BIC(pap.lca3random), NA, NA), bic3 = c(BIC(pap.lca1random2), BIC(pap.lca2random2), BIC(pap.lca3random2), NA, NA)) print(pap.bic) # plot 4 class standard plot(pap.lca4, type = "b") # plot 2 class standard plot(pap.lca2random, type = "b")# standard latent class pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) # standard latent class pap.lca3 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, cores = 1) pap.lca4 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 4, cores = 1) pap.lca5 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 5, cores = 1) # repeat with random effect with constant loading # once BIC increases fitting further models is unnecessary pap.lca1random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, random = TRUE, cores = 1) pap.lca2random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, random = TRUE, cores = 1) pap.lca3random <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, random = TRUE, cores = 1) # can also repeat fits without constant loading to give mixture of IRT models pap.lca1random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, random = TRUE, constload = FALSE, cores = 1) pap.lca2random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, random = TRUE, constload = FALSE, cores = 1) pap.lca3random2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 3, random = TRUE, constload = FALSE, cores = 1) # produce table of BIC values # shows 4 class best of standard latent class # but 2 class latent class with constant loading has better BIC pap.bic <- data.frame(bic = c(BIC(pap.lca1), BIC(pap.lca2), BIC(pap.lca3), BIC(pap.lca4), BIC(pap.lca5)), bic2 = c(BIC(pap.lca1random), BIC(pap.lca2random), BIC(pap.lca3random), NA, NA), bic3 = c(BIC(pap.lca1random2), BIC(pap.lca2random2), BIC(pap.lca3random2), NA, NA)) print(pap.bic) # plot 4 class standard plot(pap.lca4, type = "b") # plot 2 class standard plot(pap.lca2random, type = "b")
Plots the outcome probabilities for a randomLCA object, for random effects objects this can be either marginal or conditional or both. For a 2 level random effects model conditional2 will condition on the subject random effect and integrate over the period random effects. Note that plot is based on the xyplot function.
## S3 method for class 'randomLCA' plot(x, ... , graphtype = ifelse(x$random, "marginal", "conditional"), conditionalp = 0.5, classhorizontal = TRUE)## S3 method for class 'randomLCA' plot(x, ... , graphtype = ifelse(x$random, "marginal", "conditional"), conditionalp = 0.5, classhorizontal = TRUE)
x |
randomLCA object |
graphtype |
Type of graph |
conditionalp |
For a conditional graph the percentile corresponding to the random effect at which the outcome probability is to be calculated |
classhorizontal |
classes to be plotted across the page |
... |
additional parameters to xyplot |
Ken Beath [email protected]
# standard latent class with 2 classes uterinecarcinoma.lca2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, cores = 1) plot(uterinecarcinoma.lca2) uterinecarcinoma.lcarandom2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, random = TRUE, probit = TRUE, quadpoints = 61, cores = 1) # default for random effects models is marginal plot(uterinecarcinoma.lcarandom2) # default for random effects models conditional is p = 0.5 i.e. median plot(uterinecarcinoma.lcarandom2, graphtype = "conditional") # look at variability by plotting conditional probabilities at 0.05, 0.5 and 0.95 plot(uterinecarcinoma.lcarandom2, graphtype = "conditional", conditionalp = c(0.05, 0.5, 0.95))# standard latent class with 2 classes uterinecarcinoma.lca2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, cores = 1) plot(uterinecarcinoma.lca2) uterinecarcinoma.lcarandom2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, random = TRUE, probit = TRUE, quadpoints = 61, cores = 1) # default for random effects models is marginal plot(uterinecarcinoma.lcarandom2) # default for random effects models conditional is p = 0.5 i.e. median plot(uterinecarcinoma.lcarandom2, graphtype = "conditional") # look at variability by plotting conditional probabilities at 0.05, 0.5 and 0.95 plot(uterinecarcinoma.lcarandom2, graphtype = "conditional", conditionalp = c(0.05, 0.5, 0.95))
The posterior class probabililities for each observed pattern and class is determined. These are returned as a data frame together with the patterns for each observation. If class = 0 is requested then all classes are returned, otherwise only the selected class.
postClassProbs(object, class = 0)postClassProbs(object, class = 0)
object |
randomLCA object |
class |
class to be returned. Zero returns all classes. |
Extracts the corresponding data from the randomLCA object.
A data frame. The first columns of the data frame correspond to the patterns, followed by the frequency of each pattern, and then the posterior class probabilities for either the selected class or for all classes. The returned result is for the summarised data. If raw data is used, that is no frequencies, and it is required to calculated the posterior class probability for each observation then it is simply required to merge the class probabilities with the raw data, possibly removing any variable "freq" in the raw data.
Ken Beath
pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(postClassProbs(pap.lca2))pap.lca2 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 2, cores = 1) print(postClassProbs(pap.lca2))
Prints a randomLCA object. Prints summary.
## S3 method for class 'randomLCA' print(x, ...)## S3 method for class 'randomLCA' print(x, ...)
x |
randomLCA object |
... |
additional argument; currently none is used. |
Ken Beath
pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) pap.lca1 # or print(pap.lca1)pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) pap.lca1 # or print(pap.lca1)
Fit latent class models, which may include a random effect.
randomLCA(patterns, freq = NULL, nclass = 2, calcSE = TRUE, notrials = 20, random = FALSE, byclass = FALSE, quadpoints = 21, constload = TRUE, blocksize = dim(patterns)[2], level2 = FALSE, probit = FALSE, level2size = blocksize, qniterations = 5, penalty = 0.01, EMtol = 1.0e-5, verbose = FALSE, seed = as.integer(runif(1, 0, .Machine$integer.max)), cores = max(detectCores() %/% 2, 1))randomLCA(patterns, freq = NULL, nclass = 2, calcSE = TRUE, notrials = 20, random = FALSE, byclass = FALSE, quadpoints = 21, constload = TRUE, blocksize = dim(patterns)[2], level2 = FALSE, probit = FALSE, level2size = blocksize, qniterations = 5, penalty = 0.01, EMtol = 1.0e-5, verbose = FALSE, seed = as.integer(runif(1, 0, .Machine$integer.max)), cores = max(detectCores() %/% 2, 1))
patterns |
Data frame or matrix of 0 and 1 defining the outcome patterns. May also include missing values, with randomLCA using maximum likelihood to fit the models using all available data. |
freq |
Frequency for each outcome pattern, if missing this is calculated from the patterns, and the patterns are summarised to remove duplicate values. |
nclass |
Number of classes to be fitted |
calcSE |
Calculate standard errors for parameters. Useful for bootstrapping. |
notrials |
For a standard latent class model, the number of random starting values used |
random |
Include random effect? |
byclass |
Vary random effect loading(s) by class? |
quadpoints |
Number of quadrature points for adaptive quadrature |
constload |
Outcome loadings are constant for random effects model? |
blocksize |
Where a random effects (single level) model is broken into blocks, that is the loadings are repeated, this defines the size of the blocks |
probit |
Probit model for random effect? |
level2 |
Fit 2 level random effects model (further details to follow)? |
level2size |
Size of level 2 blocks if fitting 2 level models |
qniterations |
Number of Quasi-Newton iterations within each EM/adaptive cycle. Decrease if there is a failure to converge |
penalty |
penalty applied to likelihood for outcome probabilities. Shrinks outcome probabilities in slightly and can prevent extreme values. Setting penalty to 0 will produce an unpenalized fit. |
EMtol |
convergence tolerance for EM algorithm for fixed effect latent class |
verbose |
Prints fit progress if true |
seed |
Initial random seed for generating starting values. This can be set to guarantee that the fit is the same each time, including the order of the classes. |
cores |
Number of cores to be used for parallel evaluation of starting values |
The structure of the patterns is assumed to be a number of blocks of different outcomes each of level2size, allowing outcomes to be repeated. Each outcome is assumed to have it's own loading. An example is the width of the patterns is n and the level2size is n, resulting in n outcomes and therefore n loadings. Alternatively if the level2size is 1, then there are n repeats of the same outcome (but with different probabilities) with the same loading. In practice they may not be the same type of outcome, but usually will be.
The algorithm used is EM for the standard latent class and adaptive (in the sense of moving the location of the quadrature points) Gauss-Hermite quadrature for the random effects models. The number of quadrature points defaults to 21.
randomLCA object This contains
fit |
Fit object from optim |
nclass |
Number of classes |
classp |
Class probabilities |
outcomep |
Outcome probability |
lambdacoef |
Loadings |
se |
Standard errors corresponding to results returned by optim |
np |
Number of parameters |
nobs |
Number of observations in total |
logLik |
log likelihood for fitted model |
penlogLik |
Penalised log likelihood for fitted model |
observed |
Observed numbers corresponding to each pattern |
fitted |
Fitted number corresponding to each pattern |
deviance |
Deviance |
classprob |
Posterior class probability for each pattern |
bics |
BIC obtained for each trial when fitting initial latent class models |
call |
call to randomLCA |
random |
random parameter to randomLCA |
constload |
constload parameter to randomLCA |
level2 |
level2 parameter to randomLCA |
level2size |
level2size parameter to randomLCA |
byclass |
byclass parameter to randomLCA |
probit |
probit parameter to randomLCA |
quadpoints |
quadpoints parameter to randomLCA |
blocksize |
blocksize parameter to randomLCA |
freq |
frequency of each pattern |
qniterations |
qniterations parameter to randomLCA |
penalty |
penalty parameter to randomLCA |
In the returned object there a fields for patterns and frequencies. If frequencies are not supplied then the patterns and frequencies are constructed. If frequencies are supplied then zero rows are removed. When frequencies are supplied it is assumed that the data has been simplified. The returned fitted, posterior class probabilities etc, all correspond to the simplified patterns, not to the original data.
Ken Beath
Beath KJ (2017). randomLCA: An R Package for Latent Class with Random Effects Analysis. Journal of Statistical Software, 81(13), pp. 1-25. doi:10.18637/jss.v081.i13
# standard latent class with 2 classes dentistry.lca2 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, cores = 1) # random effects model with constant random effect loading dentistry.lca2random <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, random = TRUE, constload = TRUE, probit = TRUE, cores = 1) # allow loading to vary by dentist # this is the 2LCR model from Qu et al (1996) dentistry.lca2random1 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, random = TRUE, constload = FALSE, probit = TRUE, cores = 1)# standard latent class with 2 classes dentistry.lca2 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, cores = 1) # random effects model with constant random effect loading dentistry.lca2random <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, random = TRUE, constload = TRUE, probit = TRUE, cores = 1) # allow loading to vary by dentist # this is the 2LCR model from Qu et al (1996) dentistry.lca2random1 <- randomLCA(dentistry[, 1:5], freq = dentistry$freq, nclass = 2, random = TRUE, constload = FALSE, probit = TRUE, cores = 1)
Extracts the Empirical Bayes estimates of the random effects.
## S3 method for class 'randomLCA' ranef(object, ...)## S3 method for class 'randomLCA' ranef(object, ...)
object |
randomLCA object with a random effect |
... |
additional argument; currently none is used. |
A matrix with the first column containing the random effects and the second column the standard error of the random effects.
Ken Beath
pap.lca2random <- randomLCA(pap[, 1:10], freq = pap$Freq, random = TRUE, nclass = 2, cores = 1) print(ranef(pap.lca2random))pap.lca2random <- randomLCA(pap[, 1:10], freq = pap$Freq, random = TRUE, nclass = 2, cores = 1) print(ranef(pap.lca2random))
Refits an randomLCA object using new data. This is useful when fitting simulated data, for example using a bootstrap.
## S3 method for class 'randomLCA' refit(object, newpatterns, newfreq, useinit = FALSE, ...)## S3 method for class 'randomLCA' refit(object, newpatterns, newfreq, useinit = FALSE, ...)
object |
randomLCA object |
newpatterns |
the new patterns that are to be fitted using the existing model |
newfreq |
the frequencies corresponding to the patterns if required |
useinit |
use initial values from randomLCA object |
... |
additional argument; currently none is used. |
The useunit parameter determines whether the parameter estimates from the supplied model are used as initial values or whether the complete model fitting process is repeated. If the initial values are used then fitting will be faster, and the fitted classes will be similar to those in the original model. If the data was not generated from the original model there is an increased risk that the fit will not find the global maxima. For this reason when performing a bootstrap Likelihood ratio test it is better to use useinit = FALSE. However when using useinit = FALSE there may be label switching, where the estimated classes are similar, but in a different order. Unless the estimated parameters are assigned to the correct classes this will invalidate the results of a parametric bootstrap for parameter confidence intervals.
The fitted model to the new data.
Ken Beath
myocardial.lca1 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, nclass = 1, cores = 1) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1) # calculate observed lrt obslrt <- 2*(logLik(myocardial.lca2)-logLik(myocardial.lca1)) print(obslrt) nsims <- 999 # generate the simulations thesims <- simulate(myocardial.lca1, nsims) # for each simulation determin lrt simlrt <- rep(NA, nsims) for (isim in 1:nsims) { submodel <- refit(myocardial.lca1, newpatterns = thesims[[isim]]) fullmodel <- refit(myocardial.lca2, newpatterns = thesims[[isim]]) simlrt[isim] <- 2*(logLik(fullmodel)-logLik(submodel)) print(c(isim, simlrt[isim])) } # calculate p value as proportion of simulated lrt greater than observed, # corrected by adding one to numerator and denominator print((sum(simlrt >= obslrt)+1)/(nsims+1))myocardial.lca1 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, nclass = 1, cores = 1) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1) # calculate observed lrt obslrt <- 2*(logLik(myocardial.lca2)-logLik(myocardial.lca1)) print(obslrt) nsims <- 999 # generate the simulations thesims <- simulate(myocardial.lca1, nsims) # for each simulation determin lrt simlrt <- rep(NA, nsims) for (isim in 1:nsims) { submodel <- refit(myocardial.lca1, newpatterns = thesims[[isim]]) fullmodel <- refit(myocardial.lca2, newpatterns = thesims[[isim]]) simlrt[isim] <- 2*(logLik(fullmodel)-logLik(submodel)) print(c(isim, simlrt[isim])) } # calculate p value as proportion of simulated lrt greater than observed, # corrected by adding one to numerator and denominator print((sum(simlrt >= obslrt)+1)/(nsims+1))
Simulate data from a fitted randomLCA model
## S3 method for class 'randomLCA' simulate(object, nsim, seed, ...)## S3 method for class 'randomLCA' simulate(object, nsim, seed, ...)
object |
randomLCA object |
nsim |
number of data sets to be simulated |
seed |
random seed |
... |
additional optional arguments. |
Generates random data from the supplied object.
A simulated data frame or a list of simulated data frames.
Ken Beath
myocardial.lca1 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, nclass = 1, cores = 1) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1) # calculate observed lrt obslrt <- 2*(logLik(myocardial.lca2)-logLik(myocardial.lca1)) print(obslrt) nsims <- 999 # generate the simulations thesims <- simulate(myocardial.lca1, nsims) # for each simulation determin lrt simlrt <- rep(NA, nsims) for (isim in 1:nsims) { submodel <- refit(myocardial.lca1, newpatterns = thesims[[isim]]) fullmodel <- refit(myocardial.lca2, newpatterns = thesims[[isim]]) simlrt[isim] <- 2*(logLik(fullmodel)-logLik(submodel)) print(c(isim, simlrt[isim])) } # calculate p value as proportion of simulated lrt greater than observed, # corrected by adding one to numerator and denominator print((sum(simlrt >= obslrt)+1)/(nsims+1))myocardial.lca1 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, nclass = 1, cores = 1) myocardial.lca2 <- randomLCA(myocardial[, 1:4], freq = myocardial$freq, cores = 1) # calculate observed lrt obslrt <- 2*(logLik(myocardial.lca2)-logLik(myocardial.lca1)) print(obslrt) nsims <- 999 # generate the simulations thesims <- simulate(myocardial.lca1, nsims) # for each simulation determin lrt simlrt <- rep(NA, nsims) for (isim in 1:nsims) { submodel <- refit(myocardial.lca1, newpatterns = thesims[[isim]]) fullmodel <- refit(myocardial.lca2, newpatterns = thesims[[isim]]) simlrt[isim] <- 2*(logLik(fullmodel)-logLik(submodel)) print(c(isim, simlrt[isim])) } # calculate p value as proportion of simulated lrt greater than observed, # corrected by adding one to numerator and denominator print((sum(simlrt >= obslrt)+1)/(nsims+1))
Summarises the fit of a randomLCA object.
## S3 method for class 'randomLCA' summary(object, ...)## S3 method for class 'randomLCA' summary(object, ...)
object |
randomLCA object |
... |
additional argument; currently none is used. |
logLik |
Log Likelihood |
AIC |
AIC |
BIC |
BIC |
AIC3 |
AIC with penalty of 3 |
nclass |
no of classes |
probit |
link is probit |
classp |
class probabilities |
outcomep |
outcome probabilities (conditional) |
margoutcomep |
outcome probabilities (marginal), if model contains random effects |
random |
model includes random effects |
level2 |
model has 2 level hierarchy |
constload |
loadings are constant by outcome |
byclass |
lambda and tau vary by class |
lambdacoef |
lambda coefficients |
taucoef |
tau coefficients |
Ken Beath
pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) summary(pap.lca1)pap.lca1 <- randomLCA(pap[, 1:10], freq = pap$Freq, nclass = 1, cores = 1) summary(pap.lca1)
This is the data for Beath and Heller (2009).
Allergy and respiratory symptoms for infants 0 to 2 years in six month periods. Outcome is presence or absence of symptom in the six months. Original data was collected at Visits 1-7 over the 2 year period which were summarised to six month periods.
Note that these models can be slow to fit, with the "symptoms.lca2random2" model taking about 1-2 hours.
Thanks to the investigators of the CAPS study for making the data available.
symptomssymptoms
A data frame with 444 observations on the following 17 variables.
Nightcough.13Night cough in visits 1-3
Wheeze.13Wheeze in visits 1-3
Itchyrash.13Itchy rash in visits 1-3
FlexDerma.13Flexural Dermatitis in visits 1-3
Nightcough.45Night cough in visits 1-3
Wheeze.45Wheeze in visits 4-5
Itchyrash.45Itchy rash in visits 4-5
FlexDerma.45Flexural Dermatitis in visits 4-5
Nightcough.6Night cough in visit 6
Wheeze.6Wheeze in visit 6
Itchyrash.6Itchy rash in visit 6
FlexDerma.6Flexural Dermatitis in visits 1-3
Nightcough.7Night cough in visit 7
Wheeze.7Wheeze in visit 7
Itchyrash.7Itchy rash in visit 7
FlexDerma.7Flexural Dermatitis in visit 7
FreqNumber of subjects
Mihrshai et al (2001)
Mihrshahi, S., Peat, J.K., Webb, K., Tovey, R.E., Marks, G.B., Mellis, C.M. and Leeder S.R. (2001) The Childhood Asthma Prevention Study (CAPS): Design and research protocol of a randomized trial for the primary prevention of asthma. Control led Clinical Trials, 22:333–354.
Beath, K.J. and Heller, G.Z. (2009) Latent trajectory modelling of multivariate binary data. Statistical Modelling, 9(3):199–213.
symptoms.lca2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, nclass = 2, cores = 1) symptoms.lca2random <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, nclass = 2, blocksize = 4, constload = FALSE, cores = 1) symptoms.lca2random2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, level2 = TRUE, nclass = 2, level2size = 4, constload = FALSE, penalty = 0.1, cores = 1)symptoms.lca2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, nclass = 2, cores = 1) symptoms.lca2random <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, nclass = 2, blocksize = 4, constload = FALSE, cores = 1) symptoms.lca2random2 <- randomLCA(symptoms[, 1:16], freq = symptoms$Freq, random = TRUE, level2 = TRUE, nclass = 2, level2size = 4, constload = FALSE, penalty = 0.1, cores = 1)
Classification of 118 histology samples by 118 pathologists. Original classification in Holmquist et al (1967) was to one of five categories, this has been reduced to two. Analysed by a number of authors, with a random effects model in Qu et al (1996).
uterinecarcinomauterinecarcinoma
A data frame with 20 observations on the following 8 variables.
V1Pathologist 1
V2Pathologist 2
V3Pathologist 3
V4Pathologist 4
V5Pathologist 5
V6Pathologist 6
V7Pathologist 7
freqNumber of observed pattern
Qu et al (1996)
Holmquist, N.D., McMahan, C.A., and Williams, O.D. (1967) Variability in classification of carcinoma in situ of the uterine cervix. Archives of Pathology, 84, 344–345.
Qu, Y., Tan, M. and Kutner, M.H. (1996) Random effects models in latent class analysis for evaluating accuracy of diagnostic tests. Biometrics, 52, 797–810.
uterinecarcinoma.lcarandom2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, random = TRUE, probit = TRUE, quadpoints = 61, cores = 1) # LCR1 model of Que et al. This is fairly unstable and # is also slow and doesn't improve the model fit uterinecarcinoma.lcarandom2by <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, byclass = TRUE, random = TRUE, probit = TRUE, quadpoints = 71, cores = 1)uterinecarcinoma.lcarandom2 <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, random = TRUE, probit = TRUE, quadpoints = 61, cores = 1) # LCR1 model of Que et al. This is fairly unstable and # is also slow and doesn't improve the model fit uterinecarcinoma.lcarandom2by <- randomLCA(uterinecarcinoma[, 1:7], freq = uterinecarcinoma$freq, byclass = TRUE, random = TRUE, probit = TRUE, quadpoints = 71, cores = 1)