Stata The Stata listserver
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: nlcom interaction ologit

From   Michael Ingre <>
Subject   Re: st: nlcom interaction ologit
Date   Fri, 8 Apr 2005 14:00:36 +0200

On 2005-04-07, at 22:30, INAGAMI,SANAE wrote:
1)  My sample was not a random sample.  It oversampled certain
populations.  Can I still use nlcom to calculate CI's on predicted
probabilities in multilevel ordered logistic regresssions using gllamm?
If you only include only fixed effects when you calculate the CIs then it is based on an "average" subject in the population i.e. with all random effects set to zero (the group mean). If you have a strong bias in your sampling, this average subject might not be the one you think it is. Also, if you have a strong bias in your sample, the random effects might not be normally distributed which is assumed by the model. However, -gllamm- takes sampling weights that might solve your problem. See the -pweight()- option in the helpfile/manual.

You could also test if sampling predicts the estimated random effects, after the fact, by regressing group dummys or continuos variables (e.g. age if you oversampled older or younger subjects) on the latent variables. see the -geqs()- option. If they did have an effect you might have bias in your estimated variance components.

Another way is to use -gllapred u ,u- for subject specific empirical Bayes predictions of the random effects. A significant ANOVA on these predictions with your sampling groups as a factor would suggest that you have group variation that might have biased your estimates. You could also plot these estimates to assess the assumption of normality. If you use this method, be sure you only use one observation/subject.

. gllapred u , u
. egen pick = tag(subject_id)
. anova u1 sampling_group if pick
. histogram u1 if pick, normal

2)  Assuming the answer is yes...Having used the syntax that M.Ingre so
kindly provided...can I add the effect of another covariate into the
The simple answer is yes. It is possible to modify the code for ANY combination of the covariates (including random effects).

the code works...but is it ok methodology.
This is a very important question. I think it is. Some might find it odd to use a logistic model to calculate RRRs. But if I'm not mistaken, this is exactly what a multinomial logit (mlogit) does. But for other experts on this list to be able to comment, a replication of the code might help:

sysuse auto

// Create level 2 ID
gen make2 = substr(make, 1, index(make, " "))
egen pick =tag(make2)
gen id=_n if pick
sort make2 id
replace id = id[_n-1] if id == .

// Make turn>40 a dichotomous predictor variable
gen turn40 = turn > 40

// estimate gllamm model on rep78 with turn40 as predictor
gllamm rep78 turn40 , i(id) link(ologit) fam(bin) adapt

// Use local macros to temporary hold your input values
local pv1 = 1 // predictor value 1
local pv2 = 0 // predictor value 2
local cut = 11 // Choose cut point i.e. which level of the dependent variable you want to predict

// calculate log of the ratio of the two probabilities predicted from pv1/pv2
nlcom ln( (exp(`pv1' * [rep78]turn40 - [_cut`cut']_cons) / ///
(1 + exp(`pv1' * [rep78]turn40 - [_cut`cut']_cons))) / ///
(exp(`pv2' * [rep78]turn40 - [_cut`cut']_cons) / ///
(1 + exp(`pv2' * [rep78]turn40 - [_cut`cut']_cons))))

// Store the results from -nlcom-
matrix b = r(b)
matrix V = r(V)

// Estimates needs to be exponentiated to describe RRR
di "RRR: " exp(b[1,1])
di "Lower 95% CI: " exp(b[1,1] - sqrt(V[1,1])*1.96)
di "Upper 95% CI: " exp(b[1,1] + sqrt(V[1,1])*1.96)

* For searches and help try:

© Copyright 1996–2017 StataCorp LLC   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index