[R] mgcv: Plotting probabilities for binomial GAM with crossed random intercepts and factor by variable
    Jan Vanhove 
    jan.vanhove at unifr.ch
       
    Thu Jan 10 15:26:15 CET 2013
    
    
  
mgcv: Constructing probabilities for binomial GAM with crossed random 
intercepts and factor by variable
Hello,
(I'm sorry if this has been discussed elsewhere; I may not have been 
looking in the right places.)
I ran a binomial GAM in which "Correct" is modelled in terms of the 
participant's age and the modality in which the stimulus is presented 
(written vs spoken).
Participants ("Subject") and stimuli are also included as crossed random 
intercepts.
age.gam <- bam(Correct ~ Modality + s(Age, by=Modality) +
			 s(Subject, bs="re") + s(Stimulus, bs="re"),
                 data = dat, family="binomial")
Parametric coefficients:
             		Estimate Std. Error z value Pr(>|z|)
(Intercept)   		-2.242      1.121  -2.000   0.0455 *
ModalityWritten    	-1.043      1.263  -0.826   0.4087
Approximate significance of smooth terms:
                         	edf  Ref.df  Chi.sq p-value
s(Age):ModalitySpoken   	4.883   5.477   98.45  <2e-16 ***
s(Age):ModalityWritten    	5.675   6.363  126.81  <2e-16 ***
s(Subject)          		133.296 161.000 1472.35  <2e-16 ***
s(Stimulus)          		83.376  85.000 5100.59  <2e-16 ***
I would now like to plot the predicted probabilities for "Correct" == 
"yes", say for written stimuli:
plot(age.gam, select=2, shift=(-2.242-1.043), trans=plogis)
Though the overall shape of the curve is about what I expected based on 
a earlier visual exploration, the plotted probabilities are much too low 
(about 5% for values of Age where I'd expected 40%). Here are the raw 
counts:
xtabs(~ Modality + Correct, dat)
           Correct
Modality      0    1
   Spoken   4146 2693
   Written  4323 3011
Is it possible that I need to substract the values of s(Subject).1 and 
s(Stimulus).1 from coef(age.gam) to the "trans" argument, too?
Those are 0.55 and -2.86, respectively, which would provide a much 
better match between the plotted and the expected probabilities.
(But then what does the "(Intercept)" represent in this model?)
Thanks,
Jan
    
    
More information about the R-help
mailing list