Skip to content

request

9 messages · Champika Shyamalie Kariyawasam, Seth Bigelow, Steve Walker +5 more

Original
Seth Bigelow
# 
Shyam, I have found Dave Atkin's tutorial to be very helpful in a situation
that may be similar to the one you describe. I was at a loss to understand
the estimated coefficients in poisson regression before reading it. It seems
to be behind a paywall ($11.95), though it was not the last time I checked.
--Seth 

A tutorial on count regression and zero-altered count models for
longitudinal substance use data.
By Atkins, David C.; Baldwin, Scott A.; Zheng, Cheng; Gallop, Robert J.;
Neighbors, Clayton
Psychology of Addictive Behaviors, Vol 27(1), Mar 2013, 166-177.

-----Original Message-----
From: r-sig-mixed-models-bounces at r-project.org
[mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of Champika
Shyamalie Kariyawasam
Sent: Monday, May 27, 2013 11:50 PM
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] request

Hi all

I am using generalized linear mixed model fit by the laplace approximation
(family poisson) to analize my data. I have seed number (dependent variable)
as a function of study site in two lations in two countries. i ran the model
in R . But i need some assistance to interpret my data. Any source,
reference or result of previous work welcome.

thanks in advance

shyam


_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
Robert Long
# 
Hi all

I would like to obtain the design matrix for the random effects, without 
running (g)lmer first.

Could anyone help me do that ? For example, working with the sleepstudy 
dataset

sm1 <- lmer(Reaction ~ Days + (Days|Subject), data=sleepstudy)

sm1 at Zt is the transpose of the matrix I require, but I would like to 
obtain it without running lmer.

Thanks
Robert Long
Steve Walker
# 
Hi Robert,

In stable (i.e. cran) lme4, Zt can be obtained without fitting the model 
using this command:

lmer(Reaction ~ Days + (Days|Subject), data=sleepstudy, doFit = 
FALSE)$FL$trms[[1]]$Zt

In development (i.e. github) lme4, you can use this command:

lFormula(Reaction ~ Days + (Days|Subject), data=sleepstudy)$reTrms$Zt

or these two:

dv <- lmer(Reaction ~ Days + (Days|Subject), data=sleepstudy, devFunOnly 
= TRUE)
environment(dv)$pp$Zt

Cheers,
Steve
On 2013-05-28 12:44 PM, W Robert Long wrote:
Robert Long
# 
Hi Steve

That's great. Thanks a lot :-)

RL
On 28/05/2013 18:29, Steve Walker wrote:
Steven McKinney
# 
What appears to be a preprint copy is available at

   http://www.ats.ucla.edu/stat/paperexamples/atkins/Modeling_Infrequent_Counts_JFP2007_-_WEB.pdf



Steven McKinney

Statistician
Molecular Oncology and Breast Cancer Program
British Columbia Cancer Research Centre
David Winsemius
#
On May 28, 2013, at 3:49 PM, Steven McKinney wrote:

            

      
      
      
    

        
This article contains what I consider to be unfortunate advice:

"One consequence of this is that using a simple transformation allows us to interpret regression coefficients in the Poisson model as the percentage change in the expected counts:

100*[e?*? ?1] (5) where ? is the regression coefficient from the Poisson regression and ? is the units of change in

the predictor (e.g., for one unit of change in the predictor, ? = 1)."

This advice (which is repeating a misinterpretation foisted by the SAS Stats Manual) fails to recognize that the scale of interpretation is not symmetric upon "inversion" or "complementation" of the predictors. So if the predictor is a 1/0 variable, then one coding of a variable with a coefficient of log(2)  might be "interpreted" as a 100% change(say for "married"==1),  whereas the reverse coding ( alternately coded "not married" ==1)  would be "interpreted" as a 50% change. Ironically this advice occurs immediately after the distinction between linear scales and multiplicative scales. It is the multiplicative scale that invalidates that "percentage change" interpretation simply because the range from a "null"-value of 0 to the low end of possible values os "100%" whereas the range to the high end is unlimited upward.

(The correct interpretation is for the exponentiated coefficient to be seen as a relative risk or a multiplicative factor that creates a predicted count or rate relative to the Intercept or baseline estimate. There is really no need to subtract one if one realizes that a factor of 1 will not change any estimate if the result is being multiplied by a baseline value.)

If the advice were couched in more limited manner such that it were resticted to small coefficients and sufficient caveats about its approximate nature were offers, I would be less offended. It bothers me to see this further source of misinterpretation cited as an authority.

            

      
      
      
    

        

      
      
    

        
David Winsemius
Alameda, CA, USA

  
  


  


      

    

    

      
      

        


  

        

      Henrik Singmann
    

    

      
      #
    

  


  

      
Hi Robert,

A further problem could be that R is somewhat reluctant to fit models with interactions but without main effects.
I asked this on SO quite some time ago (with nice answers from regulars on this list, Ben Bolker and Joshua Wiley): http://stackoverflow.com/q/11335923/289572

Even if you specify to exclude the main effects in presence of interaction, their parameters will be there. In your data:

require(lme4)
cutsim <- read.csv(file="cutsim.csv", header=T)
cutsim$Species <- as.factor(cutsim$Species)
cutsim$Farm <- as.factor(cutsim$Farm)

# model without main effect:
model2 <- lmer(Yield ~ Species + Species:Farm + (1|Animal), data=cutsim, doFit = FALSE, REML=TRUE)
length(model2$fr$fixef)
# 500 fixed effects parameters

# model with main effect
model3 <- lmer(Yield ~ Species*Farm + (1|Animal), data=cutsim, doFit = FALSE, REML=TRUE)
length(model3$fr$fixef)
# 500 fixed effects parameters

Cheers,
Henrik


Am 28/05/2013 19:55, schrieb W Robert Long:

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

  
    

      
      
    

  


  


      

    

    
2 days later

    

      
      

        


  

        

      Ben Bolker
    

    

      
      #
    

  


  

      
Steve Walker <steve.walker at ...> writes:

            

      
      
      
    

        
It can be done even one step more directly:

mkReTrms(list(quote(Days|Subject)),sleepstudy)$Zt

the slightly convoluted form of the first argument is because
mkReTrms (see ?mkReTrms) wants a list of the random-effects
bits extracted from the formula, so it needs to be a list
of 'language' objects.  Alternatively:

mkReTrms(findbars(~(Days|Subject)),sleepstudy)$Zt

 I have a minor preference for these because they depend
a little less on the detailed form in which things are
stored in the package.

  Ben Bolker