search

dlinzer / poLCA

Polytomous Variable Latent Class Analysis (R package)
https://dlinzer.github.io/poLCA/
48 stars 17 forks source link

Bvr #3

Open daob opened 7 years ago

daob commented 7 years ago

Hi Drew,

As requested, a pull request including documentation. Function seems to work fine with covariates also.

I removed the rounding from expected counts because this can mess with downstream statistics, hope that's OK. Also I added unit tests - feel free to ignore those if you don't want them.

Best, Daniel

zh-zhang1984 commented 6 years ago

I want to make multi-group LCA with poLCA. however, I am new to the syntax and cannot figure out how to do it. based on a paper by Lanza 2013 (https://www.ncbi.nlm.nih.gov/pubmed/21318625) that is enclosed, In this paper, the author proposed two methods to investigate the effect of treatment on outcome in subgroups identified by LCA. 1) classify-analyze approach and 2) model-based approach. I want to implement the two approach with R and created following code: set.seed(8) probs <- list(matrix(c(0.6,0.2,0.2, 0.6,0.3,0.1, 0.3,0.1,0.6 ),ncol=3,byrow=TRUE), # Y1 matrix(c(0.2,0.8, 0.7,0.3, 0.3,0.7 ),ncol=2,byrow=TRUE), # Y2 matrix(c(0.3,0.6,0.1, 0.1,0.3,0.6, 0.3,0.6,0.1 ),ncol=3,byrow=TRUE), # Y3 matrix(c(0.1,0.1,0.5,0.3, 0.5,0.3,0.1,0.1, 0.3,0.1,0.1,0.5),ncol=4,byrow=TRUE), # Y4 matrix(c(0.1,0.2,0.7, 0.1,0.8,0.1, 0.8,0.1,0.1 ),ncol=3,byrow=TRUE)) # Y5 simdat <- poLCA.simdata(N=1000,probs,P=c(0.2,0.3,0.5)) trt<-as.factor(sample(c("trt","ctrl"),replace=T,size=1000)) z <- 1 - as.numeric(trt)-2simdat$trueclass+0.5as.numeric(trt)*simdat$trueclass pr <- 1/(1+exp(-z)) outcome <- rbinom(1000,1,pr) dat<-data.frame(simdat$dat,trt=trt,outcome=outcome)

classify-analyze approach

f1 <- cbind(Y1,Y2,Y3,Y4,Y5)~1 lc1 <- poLCA(f1,simdat$dat,nclass=3,nrep=5) mod<-glm(outcome~trt*as.factor(lc1$predclass), family="binomial") summary(mod)

model based approach

f2<-cbind(Y1,Y2,Y3,Y4,Y5)~outcome dat.trt<-dat[dat$trt=="trt",] dat.ctrl<-dat[dat$trt=="ctrl",] lc2.trt<-poLCA(f2,dat.trt,nclass=3,nrep=5) lc2.ctrl<-poLCA(f2,dat.ctrl,nclass=3,nrep=5) table(lc2.trt$predclass,dat.trt$outcome) prop<-rbind(ctrl=prop.table(table(lc2.ctrl$predclass,dat.ctrl$outcome),1)[4:6], trt=prop.table(table(lc2.trt$predclass,dat.trt$outcome),1)[4:6]) colnames(prop)<-c('class 1',"class 2","class 3") barplot(prop,beside =T, legend.text=c('ctrl',"trt"))

however, it appears that the model-based approach is wrong. How can I implement the model-based approach with poLCA?