[R] applying a univariate function for each response in a multivariate linear model (mlm)

[Michael Friendly]

After fitting a multivariate linear model (mlm), I'd like to be able to 
run or apply
a standard univariate stats:::*.lm function to each of the response 
variables,
within a function -- i.e., by operating on the mlm object, rather than 
re-running
the univariate models separately manually.

An example: extracting cooks.distance components (via 
stats:::cooks.distance.lm)

grain <- c(40, 17, 9, 15, 6, 12, 5, 9)        # y1
straw <- c(53, 19, 10, 29, 13, 27, 19, 30)    # y2
fertilizer <- c(24, 11, 5, 12, 7, 14, 11, 18) # x

Fertilizer <- data.frame(grain, straw, fertilizer)
# fit the mlm
mod <- lm(cbind(grain, straw) ~ fertilizer, data=Fertilizer)

# run univariate regressionsand get cooks.distance
 > (cookd.grain <- cooks.distance(lm(grain ~ fertilizer, data=Fertilizer)))
          1          2          3          4 5          6          7
3.4436e+00 4.0957e-02 2.2733e-01 4.8605e-03 1.4073e-05 2.0479e-02 
6.4192e-02
          8
4.8383e-01
 > (cookd.straw <- cooks.distance(lm(straw ~ fertilizer, data=Fertilizer)))
         1         2         3         4         5 6         7         8
2.0003953 0.0283225 0.0675803 0.1591198 0.0013352 0.0024076 0.0283225 
0.4672299

This is the result I want:

 > data.frame(cookd.grain, cookd.straw)
   cookd.grain cookd.straw
1  3.4436e+00   2.0003953
2  4.0957e-02   0.0283225
3  2.2733e-01   0.0675803
4  4.8605e-03   0.1591198
5  1.4073e-05   0.0013352
6  2.0479e-02   0.0024076
7  6.4192e-02   0.0283225
8  4.8383e-01   0.4672299

Note that if I call cooks.distance.lm directly on the mlm object, there 
is no complaint
or warning, but the result is silently WRONG:

 > # try calling cooks.distance.lm directly:  silently WRONG
 >  stats:::cooks.distance.lm(mod)
        grain      straw
1 3.4436e+00 0.51729792
2 1.5838e-01 0.02832250
3 2.2733e-01 0.01747613
4 1.8796e-02 0.15911979
5 1.4073e-05 0.00034527
6 7.9192e-02 0.00240762
7 6.4192e-02 0.00732414
8 1.8710e+00 0.46722985
 >

I realize that I can also use update() on the mlm object to re-fit the 
univariate models,
but I don't know how to extract the response names from it to do this in 
a function

 > coef(mod)  # multivariate
               grain   straw
(Intercept) -3.7524 -2.2965
fertilizer   1.4022  2.1409

 > coef(update(mod, grain ~ .))
(Intercept)  fertilizer
     -3.7524      1.4022
 > coef(update(mod, straw ~ .))
(Intercept)  fertilizer
     -2.2965      2.1409
 >


-- 
Michael Friendly     Email: friendly AT yorku DOT ca
Professor, Psychology Dept. & Chair, Quantitative Methods
York University      Voice: 416 736-2100 x66249 Fax: 416 736-5814
4700 Keele Street    Web:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA