Hello,
I'm trying to fit a special kind of proportional odds model from:
Whitehead et al. (2001). Meta-analysis of ordinal outcome using individual patient data. Statistics in medicine 20: 2243-2260. (model 2)
The data are as follows:
library(nlme) library(geepack) library(Design) library(MASS)
options(contrasts=c("contr.SAS","contr.poly"))
counts <- c(2,22,54,29,3,4,23,45,22,2,1,22,35,11,3,14,119,180,54,6,7,16,17, 10,3,13,20,24,10,1,8,24,73,52,13,21,106,175,62,17,2,13,18,7,1,3,14,19,3,0)
category <- c(rep(1:5,10))
treatment <- (rep(c(0,0,0,0,0,1,1,1,1,1,),5))
study <- gl(5,10)
percoftotal <- counts/sum(counts)
cout <- data.frame(study,treatment,category,counts,percoftotal)
The data are from five controlled clinical trials (study). The patients were given placebo (treatment=0) or treatment (treatment=1) and had a response in one of five categories.
Now I want to fit the model given by Whitehead et al. This assumes proportional odds between treatments, but stratifies by study. This means that the cut-off points associated with the distribution of the underlying latent variable for determining the response category are allowed to vary from study to study but are the same for both treatment groups within a study.
It is given by
log(Q_{ijk}/(1-Q_{ijk}) = \alpha_{ik} + \beta*x_{ij} (k=1,...,m-1)This model can be considered as arising from a latent continuos variable . Assume that the response of the j-th subject in study i is truly equal to G_{ij} although this latent reponse will never be observerd.
G_{ij} has a logistic distribution with:
Q_{ijk} = P(G_{ij} <= \alpha_{ik}) = 1/(1+exp(-(\alpha_{ik}+\beta*x_{ij}) (k=1,...,m-1)
in which Q_{ijk} is the probabilty of having a response for the j-th subject in category k or better that means p_{ij1}+....+p_{ijk}= Q_{ijk} and Q_{ijm}=1, Q_{ij1}=1.
i = 1,...,r (r=5) #number of studies
j = 1,....,n_{i} #subject j from study i
k = 1,...,m (m=5)#CategoryMy problem is how to fit this model in R although I have a SAS code from Whitehead which is availabe at http://www.rdg.ac.uk/mps/mps_home/misc/publications.htm . There it is fit by using Proc NLMIXED. The result for beta is there:
The problem is that alpha_{ik} represents the k_th intercept for the i_th study. A "normal" Proportional-Odds model
has only alpha_k which represents the k-th intercept without any associaton to the study.
I tried to fit it by using the functions of Pinheiro/ Bates nlme, nls, gnls;
polr by Venables/Ripley ;
lrm by Harrell
ordgee by Yan;
but I came to no conclusion.
I would be very happy if anyone could help me.
Thank you VERY VERY much in advance.
Kindly regards,
Henning Henke
______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
