[R] StepAIC

[Jari Oksanen]

On 29 Mar 2004, at 19:07, Christian Hennig wrote:

Dear list,

First you see the stepAIC-forward solution (fs7). The strange thing here
is that apparently not all interactions are tried for inclusion, but 
only
WQ:Lage. In particular, I think that WFL:Lage should be tried
in the last two steps, where WFL and Lage are already in the fit.
After fs7, I give the output of fs6 (backward), where all interactions 
are
tried as I have expected. (regsubsets works properly forward and
backward.)

Bob O'Hara posted a message about the same issue just one week ago. 
This seems to be a limitation in identifying names of interaction 
terms. In studying scope, A:B and B:A are regarded as different 
animals. R decides internally which way it arranges these names, and if 
the building up model has a candidate B:A, but scope has A:B, then B:A 
is not regarded as being in the scope. The net result is that step (or 
stepAIC) can be used to build interactions only with good luck.

The extreme case is that R decides to include a term (say A:B) from the 
scope, but after inclusion R decides to re-arrange its name as B:A. 
This is no longer in the scope and step ends with an error message. I 
have hoped to work out a reproducible example of this, but haven't had 
time. However, this happens with the latest devel version of vegan if 
you use methods that you shouldn¨t use (that is,  you step cca which 
you cannot do). The last step:

Step:  AIC= 125.58
  varespec ~ Al + P + K + Baresoil + P:K + P:Baresoil

               Df    AIC
+ K:Al         1 125.02
+ Zn           1 125.36
<none>           125.58
+ P:Al         1 125.60
+ Al:Baresoil  1 125.82
+ Humdepth     1 125.83
+ Mo           1 125.94
+ Mg           1 125.96
+ Mn           1 126.11
+ S            1 126.36
- P:Baresoil   1 126.43
+ N            1 126.52
+ Fe           1 126.72
+ K:Baresoil   1 126.80
+ pH           1 126.89
+ Ca           1 127.01
- P:K          1 127.14
- Al           1 128.01

Step:  AIC= 125.02
  varespec ~ Al + P + K + Baresoil + P:K + P:Baresoil + Al:K

Error in factor.scope(ffac, list(add = fadd, drop = fdrop)) :
	upper scope does not include model

Interpretation: K:Al was in the scope and it was included. However, 
after inclusion it changed into Al:K which is not in the scope, so that 
R was able to produce a model where "upper scope does not include 
model".

The secret is in the C routines which decide how to order terms in 
formulae.

cheers, jari oksanen (Oulu)