R help - May 2012 - mgcv: How to calculate a confidence interval of a ratio

Dear R-Users,

Dr. Wood replied to a similar topic before where confidence intervals were
for a ratio of two treatments (
https://stat.ethz.ch/pipermail/r-help/2011-June/282190.html). But my
question is more complicated than that one. In my case, log(E(y)) = s(x)
where y is a smooth function of x. What I want is the confidence interval
of a ratio of log[(E(y2))/E(y1)] given two fixed x values of interest. This
is complicated than two treatments because they can be modeled as a binary
variable as Dr. Wood pointed out. I am wondering if mgcv has some embedded
functions to calculate this quickly.

Best regards,
Gevin

2012-05-23

	[[alternative HTML version deleted]]

Given that this is just s(x2) - s(x1) then you can get the CI using the
type= "lpmatrix" with predict.gam. Here's an example...

library(mgcv)
## simulate some data
dat <- gamSim(1,n=400,dist="poisson",scale=.25)
## fit log-linear model...
b <- gam(y~s(x0)+s(x1)+s(x2)+s(x3),family=poisson,
          data=dat,method="REML")
## data at which predictions to be compared...
pd <- data.frame(x0=c(.2,.3),x1=c(.5,.5),x2=c(.5,.5),
                  x3=c(.5,.5))
## log(E(y_1)/E(y_2)) = s(x_1) - s(x_2)
Xp <- predict(b,newdata=pd,type="lpmatrix")
## ... Xp%*%coef(b) gives log(E(y_1)) and log(E(y_2)),
## so the required difference is computed as...
diff <- (Xp[1,]-Xp[2,])
dly <- t(diff)%*%coef(b) ## required log ratio (diff of logs)
se.dly <- sqrt(t(diff)%*%vcov(b)%*%diff) ## corresponding s.e.

dly + c(-2,2)*se.dly ## 95%CI




On 23/05/12 15:37, Gevin Brown wrote:
> Dear R-Users,
>
> Dr. Wood replied to a similar topic before where confidence intervals were
> for a ratio of two treatments (
> https://stat.ethz.ch/pipermail/r-help/2011-June/282190.html). But my
> question is more complicated than that one. In my case, log(E(y)) = s(x)
> where y is a smooth function of x. What I want is the confidence interval
> of a ratio of log[(E(y2))/E(y1)] given two fixed x values of interest. This
> is complicated than two treatments because they can be modeled as a binary
> variable as Dr. Wood pointed out. I am wondering if mgcv has some embedded
> functions to calculate this quickly.
>
> Best regards,
> Gevin
>
> 2012-05-23
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Simon Wood, Mathematical Science, University of Bath BA2 7AY UK
+44 (0)1225 386603               http://people.bath.ac.uk/sw283