[R] fit a nonlinear model using nlm()
I am trying to fit a nonlinear model using nlm().
My application of nlm() is a bit complicated.
Here is the story behind the model being fit:
The observer is trying to detect a signal corrupted by noise.
On each trial, the observer gets stim=signal+rnorm().
In the simulation below I have 500 trials. Each row of stim is a new trial.
On each trial, if the cross-correlation between the stim and
the signal is above some criterion level (crit=.5 here), the
observer says signal (resp=1), else he says no signal
(resp=0).
So the situation is this: I know resp and stim for each trial.
I would like to estimate crit and signal from these.
(You might say that I already know the signal and crit. But the
idea here is that the observer cross-correlates the stim with an
internal template that will not be identical to the actual signal.
I want to estimate this template. Also the observer's crit may
differ from the correct one.)
In the code below, please help me get the f() and nlm() bits right.
I want to estimate signal and crit given stim and resp.
Thanks very much for any help!
Bill
x-1:100
con-.1
signal-con*cos(2*pi*3*x/length(x))
crit-.5
noisesd-.1
# each row is a new stim (trial). 500 trials
resp-array(dim=500)
stim-matrix(nrow=500,ncol=100)
for (i in 1:500)
{
stim[i,]-signal+rnorm(n=length(signal), mean=0, sd=noisesd)
if (sum(stim[i,]*signal)crit) (resp[i]-1) else (resp[i]-0)
}
f-function(signalest)
{
r-array(dim=500)
for (i in 1:500)
{
r[i]-sum(stim[i,]*signalest)critest
}
sum((r-resp)^2)
}
fit-nlm(f, stim[1,])
__
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
and provide commented, minimal, self-contained, reproducible code.
[R] dose-response on a grid
I have the following problem. I have measured a dose response curve (binary response, continuous dose) on a grid of x,y positions. I would like to produce a grey-level plot that shows the LD50 at each (x,y) position. I am thinking that I have to do something like fit-glm(resp ~ x*y + dose, family = binomial) Corrections welcome. But from here I don't know how to get LD50, and certainly not at each x,y, position. Thanks very much for any help. Bill __ 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 and provide commented, minimal, self-contained, reproducible code.
[R] lme() F-values disagree with aov()
I have used lme() on data from a between-within subjects experiment. The correct ANOVA table is known because this is a textbook example (Experimental Design by Roger Kirk Chapter 12: Split-Plot Factorial Design). The lme() F-values differ from the known results. Please help me understand why. d-read.table(kirkspf2.dat,header=TRUE) for(j in 1:4) d[,j] - factor(d[,j]) ### Make vars into type factor ##lme() results library(nlme) fit-lme(y~a*b*c,random=~1|s, data=d) anova(fit) ##correct anova table ##subjects are nested within a; a between, b c within fit2-aov(y ~ a*b*c + Error(s/(c*b)), data=d) summary(fit2) I suspect I need a different random=... statement in lme(). Thanks very much for any help Bill The data file is attached -- kirkspf2.dat Here it is again: s a c b y 1 1 1 1 3 1 1 1 2 7 1 1 2 1 4 1 1 2 2 7 2 1 1 1 6 2 1 1 2 8 2 1 2 1 5 2 1 2 2 8 3 1 1 1 3 3 1 1 2 7 3 1 2 1 4 3 1 2 2 9 4 1 1 1 3 4 1 1 2 6 4 1 2 1 3 4 1 2 2 8 5 2 1 1 1 5 2 1 2 5 5 2 2 1 2 5 2 2 2 10 6 2 1 1 2 6 2 1 2 6 6 2 2 1 3 6 2 2 2 10 7 2 1 1 2 7 2 1 2 5 7 2 2 1 4 7 2 2 2 9 8 2 1 1 2 8 2 1 2 6 8 2 2 1 3 8 2 2 2 11 __ 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme() F-values disagree with aov()
Your lme statement is OK. To get the usual split-plot anova, your aov statement should be fit2 - aov(y ~ a*b*c + Error(s), data = d) No, this gives wrong F-values. By wrong I mean it does not agree with the published table. Table 12.10-2, page 559: Number of obs = 32 R-squared = 0.9920 Root MSE = .559017 Adj R-squared = 0.9589 Source | Partial SSdf MS F Prob F ---+ Model | 233.62525 9.345 29.90 0.0002 | a | 3.125 1 3.125 2.00 0.2070 s|a | 9.375 6 1.5625 ---+ b | 162.00 1 162.00 199.38 0. a*b | 6.125 1 6.125 7.54 0.0335 b*s|a | 4.875 6 .8125 ---+ c | 24.50 1 24.50 61.89 0.0002 a*c | 10.125 1 10.125 25.58 0.0023 c*s|a | 2.375 6 .39583 ---+ b*c |8.00 18.00 25.60 0.0023 a*b*c | 3.125 1 3.125 10.00 0.0195 | Residual | 1.875 6 .3125 ---+ Total | 235.5031 7.59677419 Bill __ 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme() F-values disagree with aov()
Thanks very much Peter! Your lme statement is OK. To get the usual split-plot anova, your aov statement should be fit2 - aov(y ~ a*b*c + Error(s), data = d) No, this gives wrong F-values. By wrong I mean it does not agree with the published table. Well, it's the model that is equivalent to your lme() model Yes, I thought my lme() model was wrong but couldn't figure out how to do it properly. Thing is that you want to add random effects of s:b and s:c, which are crossed factors, so somewhat tricky to code with lme() (this sort of thing is easier in lmer() from the lme4 packages). I am happy to do it that way if you show me how... The generic way to handle this in lme() is via something like random=list(s=pdBlocked(list( pdIdent(~1), pdIdent(~b-1), pdIdent(~c-1 I see. Thanks very much Peter!! Bill __ 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] between-within anova: aov and lme
Well nobody answered :-( Nobody on R-help is doing anovas I guess -- I don't blame them! (It's just for aggies.) In the absence of any response and for no good reason I am doing: fitn1 - aov(amplitude ~ stereo*site*stimulus + Error(subject), stereon1) This is Bill Venables's way. And when the data are unbalanced I am doing: lme(amplitude ~ site+stimulus+stereo*stimulus, random=~1|subject, method=ML, stereon1) And I have no clue why. Every discussion of between-within ANOVA I have read (practical or mathematical) is either vacuous or opaque... Cheers Bill I have 2 questions on ANOVA with 1 between subjects factor and 2 within factors. 1. I am confused on how to do the analysis with aov because I have seen two examples on the web with different solutions. a) Jon Baron (http://www.psych.upenn.edu/~baron/rpsych/rpsych.html) does 6.8.5 Example 5: Stevens pp. 468 - 474 (one between, two within) between: gp within: drug, dose aov(effect ~ gp * drug * dose + Error(subj/(dose*drug)), data=Ela.uni) b) Bill Venables answered a question on R help as follows. - factor A between subjects - factors B*C within subjects. aov(response ~ A*B*C + Error(subject), Kirk) An alternative formula would be response ~ A/(B*C) + Error(subject), which would only change things by grouping together some of the sums of squares. --- SO: which should I do aov(response ~ A*B*C + Error(subject), Kirk) aov(response ~ A/(B*C) + Error(subject), Kirk) aov(response ~ A*B*C + Error(subject/(B*C)), Kirk) 2. How would I do the analysis in lme()? Something like lme(response~A*B*C,random=~1|subject/(B*C))??? Thanks very much for any help! Bill Simpson __ 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 and provide commented, minimal, self-contained, reproducible code. __ 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 and provide commented, minimal, self-contained, reproducible code.
[R] between-within anova: aov and lme
I have 2 questions on ANOVA with 1 between subjects factor and 2 within factors. 1. I am confused on how to do the analysis with aov because I have seen two examples on the web with different solutions. a) Jon Baron (http://www.psych.upenn.edu/~baron/rpsych/rpsych.html) does 6.8.5 Example 5: Stevens pp. 468 - 474 (one between, two within) between: gp within: drug, dose aov(effect ~ gp * drug * dose + Error(subj/(dose*drug)), data=Ela.uni) b) Bill Venables answered a question on R help as follows. - factor A between subjects - factors B*C within subjects. aov(response ~ A*B*C + Error(subject), Kirk) An alternative formula would be response ~ A/(B*C) + Error(subject), which would only change things by grouping together some of the sums of squares. --- SO: which should I do aov(response ~ A*B*C + Error(subject), Kirk) aov(response ~ A/(B*C) + Error(subject), Kirk) aov(response ~ A*B*C + Error(subject/(B*C)), Kirk) 2. How would I do the analysis in lme()? Something like lme(response~A*B*C,random=~1|subject/(B*C))??? Thanks very much for any help! Bill Simpson __ 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 and provide commented, minimal, self-contained, reproducible code.
