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

From: Martin Henry H. Stevens <HStevens_at_muohio.edu>
Date: Fri 16 Jun 2006 - 02:05:40 EST

Hi John,
I think a solution is to
1. recode A and B as a single factor, AB, with four levels, 2. define each fixed effect as a function of AB minus the intercept (e.g. ed50 ~ as.factor(AB)-1).
3. extract the tTable as a data.frame with summary(model)$tTable.

I will be interested to see what other folks suggest. Cheers,
Hank

  and then run the Probably not the best, nor the worst solution, might be to recode A and B
On Jun 15, 2006, at 11:36 AM, John S. Walker wrote:

> 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.
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-
> guide.html

Dr. M. Hank H. Stevens, Assistant Professor 338 Pearson Hall
Botany Department
Miami University
Oxford, OH 45056

Office: (513) 529-4206
Lab: (513) 529-4262
FAX: (513) 529-4243

http://www.cas.muohio.edu/~stevenmh/
http://www.muohio.edu/ecology/
http://www.muohio.edu/botany/

"E Pluribus Unum"

R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Fri Jun 16 02:24:17 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Fri 16 Jun 2006 - 04:12:58 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.