[R] Sum of the deviance explained by each term in a gam model does not equal to the deviance explained by the full model.

Simon Wood s.wood at bath.ac.uk
Fri Nov 11 15:10:52 CET 2011


It's the same problem as with any regression model where the predictors 
are not strictly orthogonal. You can make the explained deviances per 
term sum to the whole model explained deviance by dropping terms 
sequentially, but then the explained deviance per term depends on the 
order in which you drop terms. The alternative in the code you pasted in 
below is at least independent of the ordering.

I suppose you could make things add up in the way you want by 
symmetrising the computation like this...

## reduce model one way....
dev.1<- (deviance(b1)-deviance(b))/deviance(b0) ## prop expl by s(x2)
dev.2<- (deviance(b0)-deviance(b1))/deviance(b0) ## prop expl by  s(x1)

## reduce it the other way...
dev.2[2]<- (deviance(b2)-deviance(b))/deviance(b0) ## prop expl by s(x2)
dev.1[2]<- (deviance(b0)-deviance(b2))/deviance(b0) # prop expl by s(x1)

## average the alternatives
dev.1 <- mean(dev.1)
dev.2 <- mean(dev.2)

## so now explained deviance adds up...
dev.1+dev.2
summary(b)$dev.expl

... but I'm not convinced that this really adds anything useful, and it 
gets very tedious as the number of model terms increases...

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


On 10/11/11 14:02, Huidong TIAN wrote:
> Dear R users,
>      I read your methods of extracting the variance explained by each
> predictor in different places. My question is: using the method you
> suggested, the sum of the deviance explained by all terms is not equal to
> the deviance explained by the full model. Could you tell me what caused
> such problem?
>
>>   set.seed(0)
>>   n<-400
>>   x1<- runif(n, 0, 1)
>>   ## to see problem with not fixing smoothing parameters
>>   ## remove the `##' from the next line, and the `sp'
>>   ## arguments from the `gam' calls generating b1 and b2.
>>   x2<- runif(n, 0, 1) ## *.1 + x1
>>   f1<- function(x) exp(2 * x)
>>   f2<- function(x) 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10
>>   f<- f1(x1) + f2(x2)
>>   e<- rnorm(n, 0, 2)
>>   y<- f + e
>>   ## fit full and reduced models...
>>   b<- gam(y~s(x1)+s(x2))
>>   b1<- gam(y~s(x1),sp=b$sp[1])
>>   b2<- gam(y~s(x2),sp=b$sp[2])
>>   b0<- gam(y~1)
>>   ## calculate proportions deviance explained...
>>   dev.1<- (deviance(b1)-deviance(b))/deviance(b0) ## prop explained by
> s(x2)
>>   dev.2<- (deviance(b2)-deviance(b))/deviance(b0) ## prop explained by
> s(x1)
>>
>>   dev.1 + dev.2
> [1] 0.6974949
>>   summary(b)$dev.expl
> [1] 0.7298136
>
> I checked the two models (b1&  b2), found the model coefficients are
> different with model b, so I feel it could be the problem.
>
> wish to hear your comments.
>
> Huidong Tian
>
> 	[[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.
>



More information about the R-help mailing list