
First of all, thank you for replying me.

I believe that using a statistical software (like R) and understanding the statistical issues beyond the software are two concepts with a strong link, but I understand that your scope is providing information on the way R works (so how to use it). For such a reason I know my questions were not really pertinent (or not at all!).
The truth is that it's summer-time and I have serious problems in getting some statistical advice that I really need. I'm studying and reading a lot and in fact I'm improving my poor statistical knowledge every day (I'm a biologist) but, you know, it's not easy...    
Your response has been very useful for me and I appreciate it a lot. Below I give some commentary (to better explain things and to thank you again) and I make some very short question.
I appreciate your effort although you don't answer me again. I hope not to seem hypocrite (I'm not) but I'm very grateful to you.
I promise to write again only for more pertinent questions.

Best wishes


> Date: Fri, 10 Jul 2009 17:45:50 -0700
> From: bolker@ufl.edu
> To: r-help@r-project.org
> Subject: Re: [R] generalized linear model (glm) and "stepAIC"
> 
> 
> 
>   If you possibly can, you should get some local statistical advice.
> There are a number of pitfalls involved in what you're doing.
> Frank Harrell's book is a good reference, but may be too advanced
> if you are a beginner.  You are right on the edge of having too little
> data for what you want to do (the rule of thumb is that you should
> have at least 10-20 responses per predictor), and stepwise regression is
> known
> to inflate p-values (see e.g. Whittingham et al J. Animal Ecology 2006).
> That said, here are some specific responses:
> 
> 
> Simone Santoro wrote:
> > 
> > 
> > I have 12 response variables (species growth rates) and two
> > environmental factors that I want to test to find out a possible
> > relation.
> > 
> > The sample size is quite small: (7<n<12, depending on each species-case).
> > 
> > I performed a Shapiro test (shapiro.test) to test for normal
> > distribution of the responses variables and they were normally
> > distribuited 10 times (over 12 possible, i.e. 12 response variables).
> > 
> 
>   The Shapiro test is probably not very powerful for such a small
> data set -- i.e., the data could be non-normal (in fact it almost
> certainly *is* non-normal) but the deviation
> is not detectable ...  where do your growth rates come from?
> Can you make a guess at their probable distribution?
>
> 
>>> The growth rates are calculated as  ÄXt, where ÄXt = (Xt + 1) - Xt , Xt is loge (Nt), and Nt is the population size at time t. 
>>> I use it and not directly population size because I found in a few cases (species population size trend) the existence of autocorrelation (time lag= >>> 1), nevertheless the "ÄXt"didn't >>> show autocorrelation and was equivalent to my purpose: investigating if "x1" or "x2" affected to the population dynamic of these species. I would expect that "ÄXt" would be normally >>> distribuited. 
>
> 
> > I performed a Generalized Linear Model in R-software (MASS package),
> > and I selected models by automatic backward stepwise (stepAIC
> > procedure) considering as the starting model the one with the additive
> > effects of both the factors. This is the case for six species growth
> > rates (six cases) but for the others six I tested the effect of just
> > one factor ("x2", see below) using just the "glm" procedure.
> > 
> 
> If you're assuming normality you don't need glm(), just lm(),
> although glm() will give you the same answer less efficiently
> (a common confusion, probably due in part to SAS's PROC GLM,
> which is about the same as lm()).
>
>
>>> Ok, I understand.
>
>
> Why different procedures for different cases?
>
> 
>>> I don't understand if you are suggesting to me to use different procedures for different cases or if you are asking me why I used different procedures for different cases, in such case: I >>> didn't. 
>
> 
> > So, my object containing the data is called "data" and, this is the editor
> > for the first species (sp1):
> > 
> > GLM1<-glm(growth.sp1~x1+x2,family=gaussian, data)
> > 
> > MOD.SELECTION<-stepAIC(GLM1, trace=TRUE) 
> > 
> > summary(MOD.SELECTION)
> > 
> > 
> > Here I attach an example of one of these analyses and after I finally
> > give you my questions (I hope not to be too long-winded!!):
> > 
> > 
> > 
> >> sp1.starting.model<-glm(sp1~x1+x2,family=gaussian, data)
> > 
> >> sp1.back<-stepAIC(sp1.starting.model, trace=TRUE)
> > 
> > Start:  AIC=63.6
> > 
> > sp1 ~ x1 + x2
> > 
> >        Df Deviance     AIC
> > 
> > - x2     1   73.490  61.801
> > 
> > <none>      72.278  63.602
> > 
> > - x1    1  122.659  67.949
> > 
> > 
> > 
> > Step:  AIC=61.8
> > 
> > gpf ~ x1
> > 
> > 
> > 
> >        Df Deviance     AIC
> > 
> > <none>      73.490  61.801
> > 
> > - x1    1  126.400  66.309
> > 
> >> summary(sp1.back)
> > 
> > 
> > 
> > Call:
> > 
> > glm(formula = sp1 ~ x1, family = gaussian, data = data)
> > 
> > 
> > 
> > Deviance Residuals: 
> > 
> >      Min        1Q    Median        3Q       Max  
> > 
> > -5.04833  -1.15233  -0.06802   0.81325   5.11464  
> > 
> > 
> > 
> > Coefficients:
> > 
> >             Estimate Std. Error t value Pr(>|t|)  
> > 
> > (Intercept) -7.62399    3.11127  -2.450   0.0342 *
> > 
> > x1           0.20595    0.07675   2.683   0.0230 *
> > 
> > ---
> > 
> > Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 
> > 
> > 
> > 
> > (Dispersion parameter for gaussian family taken to be 7.348965)
> > 
> > 
> > 
> >     Null deviance: 126.40  on 11  degrees of freedom
> > 
> > Residual deviance:  73.49  on 10  degrees of freedom
> > 
> >   (1 observation deleted due to missingness)
> > 
> > AIC: 61.801
> > 
> > 
> > 
> > Number of Fisher Scoring iterations: 2
> > 
> > 
> > 
> 
>  You would probably be better off just doing summary() and
> looking at the p-values of the two predictors (if you must ...)
> 
> Why are you using AIC if you're interested in testing
> relationships rather than prediction?
>
> 
>>> So, by reading the Whittingham et al. paper (thank you very much) and reading your commentaries I undertand I would be better off using the "full" (just two predictors) model and >>> taking in account the p-values of such a model (not using the stepAIC procedure), isn't it?  
>
>  
> > THE QUESTIONS:
> > 
> > 1) Can I trust in the existence of such statistical relation? I mean: is
> > there a way to know the power of this test in R?
> > 
> 
> There are power tests in R, but I don't know if there are any specifically
> for this case (two-predictor regression).  Remember that power applies to
> the probability of type II (falsely failing to reject null hypothesis)
> errors.
>
>
>>> Ok, on the other hand, I suppose that the small sample size makes the existence of a statistical relation between the predictor and the response variable even more reliable, isn't it?
> 
> 
> > 2) I decided to use always "family=gaussian" because I have also
> > negative values in my response variable and I cannot manage it in a
> > different way. In fact I was not able to use as link function a
> > "negative binomial" as I previously did in SAS because of negative
> > values of response variable (as R "told" me when I tried)
> > 
> 
> Is this a question?  As above, glm() with gaussian family and
> identity (default) link is equivalent to lm().
> 
>
>>> Yes, now I understand and I shame because I'm aware it is a very basic statistical issue (I'm sorry!). But, if I strongly believe the response variable is normally distribuited, although >>> the small sample size makes difficult to test its normality, can I use lm() without testing for normality? In other words: can I trust on logical basis that the statistical population beyond >>> the sample would be normally distribuited and consequently using lm()? 
>
> 
> > 3) How should I interpret the dispersion value R give me (in the case
> > reported it was "7.3")? I mean, what range of values (if it does exist)
> > I would expect to make the result reliable in the case of
> > "family=gaussian" (I'm not interested in prediction but just in finding
> > a statistical relation)?
> > 
> 
> In the case of the Gaussian, the dispersion is just the variance.
> If you're worried about overdispersion, that's really just a problem
> with models (binomial, Poisson) that assume a fixed relationship
> between mean and variance.
> 
>
>>> I was looking a lot for this information but I didn't find it: Thank you very much!!
>
>
>   If you really were going to use AIC, you should strongly
> consider AICc -- corrected for small sample sizes ...
> 
>   Ben Bolker
> 
> -- 
> View this message in context: http://www.nabble.com/generalized-linear-model-%28glm%29-and-%22stepAIC%22-tp24433331p24436223.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> ______________________________________________
> R-help@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.

_________________________________________________________________


	[[alternative HTML version deleted]]


