[R] Repost: Estimation when interaction is present: How do I get get the parameters from nlme?

John S. Walker jsw9c at uic.edu
Thu Jun 15 17:36:39 CEST 2006


Gday,

This is a repost since I only had one direct reply and I remain 
mystified- This
may be stupidity on my part but it may not be so simple.




In brief, my problem is I'm not sure how to extract parameter 
values/effect sizes from a nonlinear
regression model with a significant interaction term.

My data sets are dose response curves (force and dose) for muscle that 
also have two treatments applied
Treatment A (A- or A+) and Treatment B (B-/B+). A single muscle was 
used for each experiment - a full dose response curve and one treatment 
from the matrix A*B (A-/B-, A+/B-, A-/B+ and A+,B+). There are 8 
replicates for each combination of treatments
We fit a dose response curve to each experiment with parameters upper, 
ed50 and slope; we expect treatment A to change upper and ed50. We want 
to know if treatment B blocks the effect of treatment A and if so to 
what degree.
This is similar to the Ludbrook example in Venables and Ripley, however 
they only had one treatment and I have two.

my approach

The dataframe is structured like this:
expt 	treatA 	treatB 	dose 	force.
1	-		-		0.1		20
1	-		-		0.2		40
...
4	+		+		0.1		20
4	+

I used a groupedData object:  mydata=groupedData(force ~ dose | expt)

I used an nlme obect to model the data as follows (pseudocode):

myfit.nlme <- nlme(force ~ ss_tpl(dose, upper, ed50,slope), 
fixed=list(ed50~factor(treatA)*factor(treatB)))


The function ss_tpl  is a properly debugged and fully functional
selfstarting three parameter logistic function that I wrote- no problem 
here. In my analysis
I also included fixed terms for the other fit parameters; upper and 
slope, but my main problem  is with  the
ed50 so that's all I've included here.

Running an anova on the resulting object (anova(myfit.nlme) I found the 
A -/B- (control) to
be significantly different from zero, treatment A was significantly 
different, treatment B had no significant
effect  and there was a significant interaction between treatment A and 
treatment B.

The interaction term is likely to be real. The treatments are on
sequential steps in a pathway and treatment B may be blocking the
effect of treatment A, i.e. treatment B alone has no effect because it
blocks a pathway that is not active, treatment A reduces force via this
pathway and treament B therefore blocks the effect of treatment A when
used together.

 From what I understand, please correct me if I'm wrong, the parameter 
estimates from summary(model.nlme) are not correct for main effects if 
a significant interaction is present. For example in my data treatment 
B alone has no signifcant effect in the anova but the interaction term 
A:B is significant. I believe The summary estimate for B is the 
estimate across all levels of A. What I want to do is pull out the 
estimate for B when A is not present. I suppose I can do it manually 
from the list of coefficients from nls or fit a oneway model with 
treatment levels A, B, AB. But I was kind of hoping there was some 
extractor function.

The reason I need this is that the co-authors want to include a table 
of  parameter values  with std errs or confidence intervals ala:

Treat	upper	ed50	slope
A-/B-		x		x		x		<- shows value for comparison to control studies
A+/B-	x		x		x		<-Shows A is working0
A-/B+	x		x		x		<- Shows B has no effect alone
A+/B	+	x		x		x 		<-shows B blocks A (not necessarily total)
	
So back to my question,How do I extract estimates of the parameters 
from my model object for a
specific combination of factors including the interaction term.
   i.e. what is the ed50 (and std err) for A-/B-, A+/B-, A-/B+, A+/B+ ?


I think this is a fair question and one that many biomedical scientists 
would need.



More information about the R-help mailing list