[R] Calculate Total CORRECTED SS for non linear regression

[Michael Eisenring]

Hello

I started to work on nonlinear regressions recently,

I want to compare a set of models and would like to use an equivalent of R2
for  nonlinear regression models.

I learned that a good alternative is to calculate the following corrected
formula:

 

 1 - (Residual SS/Total Corrected SS)

 

This provides a measure of how well the model fits versus just the overall
mean for the response variable.

 

With some help I could program the formula for my model and , but I am not
quite sure if it is correct and what exactly the total CORRECTED SS is
(compared to the Total SS)

Can anyone tell me how to calculate the Total Corrected SS in R and how it
can be implemented in my code?

 

CODE and INPUT:

 

 


# REGRESSION NON-LINEAR REGRESSION MODEL
FIT------------------------------------
dta<-read.csv("input.csv",header=T, sep = ",")
dput(dta)
head(dta)

 

# loading packages: analysis of mixed effect models
library(nls2)#model

 

 


#EQUATION 1: y~yo+a*(1-b^x)
#Aim: fit equation to data: y~yo+a*(1-b^x)
# y =Compound (from my data set) x= Damage (from my data set)
#The other 3 parameters are unknown: yo=Intercept, a= assymptote ans b=slope

plot(Compound~Damage, dta)
# Looking at the plot, 0 is a plausible estimate for y0:
# a+y0 is the asymptote, so estimate about 4000;
# b is between 0 and 1, so estimate .5
dta.nls <- nls(Compound~y0+a*(1-b^Damage), dta,
               start=list(y0=0, a=4000, b=.5))

xval <- seq(0, 10, 0.1)
lines(xval, predict(dta.nls, data.frame(Damage=xval)))
profile(dta.nls, alpha= .05)


summary(dta.nls)


#Calculation R2 Variation expl: total Variation--> Is deviance(dta.nls) the
same as Total Corrected SS???
CompoundSS <- sum((dta$Compound - mean(dta$Compound))^2)
R2 <- deviance(dta.nls)/CompoundSS
R2
#R2=0.5723023

 

 

 

structure(list(Treatment = structure(c(1L, 2L, 3L, 1L, 2L, 3L, 
4L, 1L, 2L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 
1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 
1L, 2L, 3L, 4L, 1L, 2L, 4L), .Label = c("1_2d", "3_2d", "9_2d", 
"C"), class = "factor"), Compound = c(1036.331811, 4171.427741, 
6039.995102, 4140.242559, 4854.985845, 6982.035521, 948.2418407, 
3618.448997, 3130.376482, 1180.171957, 1500.863038, 4576.787021, 
5629.979049, 3589.187889, 2508.417927, 1989.576826, 5972.926124, 
450.7205451, 1120.955, 3470.09352, 3575.043632, 349.0864019, 
1013.807628, 910.8879471, 3743.331903, 592.3403778, 1517.045807, 
1504.491931, 3736.144027, 723.885643, 1782.864308, 1414.161257, 
3723.629772, 2005.919344, 4198.569251, 2228.522959, 3322.115942, 
720.9785449, 2874.651764, 2287.228752, 5654.858696, 1247.806111, 
2547.326207, 2608.716056, 1079.846532), Damage = c(0.4955, 1.516, 
4.409, 0.491, 2.3035, 3.51, 0, 0.4435, 1.573, 0, 0.142, 2.171, 
4.023, 0, 0.6925, 1.989, 5.683, 0, 0.756, 2.129, 9.437, 0, 0.578, 
2.966, 4.7245, 0, 1.0475, 1.62, 5.568, 0, 0.8295, 2.411, 7.272, 
0, 0.4035, 2.974, 8.043, 0, 0.6965, 1.313, 5.681, 0, 0.5895, 
2.559, 0)), .Names = c("Treatment", "Compound", "Damage"), class =
"data.frame", row.names = c(NA, 
-45L))

 

 

 


	[[alternative HTML version deleted]]