[R] Lattice: raw data and prediction of a non linear fitted function
François Collin
fanch.collin at gmail.com
Thu Aug 6 12:34:10 CEST 2015
Hi,
Thank you for your answer. However, I think you misunderstood my problem.
The data and fit I provided are just here to illustrate 2 points:
- I have few raw data
- and a fitted model I want to represent.
It is a simulated dataset, no more, which match to my requirement: a non
linear response, which 2 unknown parameters b2 and b3 (?SSgompertz) have
values according to their group (char = { A, B }). The fit I proposed take
this into account as my formula is:
y ~ SSgompertz(time, Asym = 100, b2[char], b3[char]).
So I do not think two equations are relevant (or is out of the bounds of my
problem). But again, the quality of the fit is not my concern up to this
point, these are just required to illustrate my problem.
My question is about graphical representation using lattice. I look for the
efficient way to represent: my raw data and a curve which represents my
fitted function for two groups. But, as I want a curve, I have to increase
the number of points for the prediction to obtain a smooth representation.
The second part of your answer is more relevant for me:
================
xyplot(y ~ time, data = foo,
groups = char,
panel = panel.superpose,
panel.groups = function(x, y, ...,group.number){
panel.xyplot(x,y, ...);
xy.nls <-
nls(y ~ SSgompertz(x, Asym = 100, b2, b3),
start = list( b2 =7, b3 = 0.9))
#print(summary(xy.nls))
if (group.number == 1){
xy.pre <- predict(xy.nls, newdata = list(time =
seq(0,55,5)))
print(xy.pre)
llines( seq(5,60,5), xy.pre, col = "blue")
} else {
xy.pre <- predict(xy.nls, newdata = list(time =
seq(5,60,5)))
print(xy.pre)
llines( seq(5,60,5), xy.pre, col = "magenta")
}
})
==============
However, I tried it, and the results is very weird: instead of having a
smoother representation, I get a repetition of the initial pattern which
looks periodic, and draw according the number of points of the predict but
not using the predicted y.
Thus my problem is still unsolved:
- How can I use the group argument inside xyplot to make new predictions
inside panel
argument ?
- How can I use this prediction inside the panel argument to draw the fit
per group ?
- My constraints: on the same graph I want my raw data and my predictions
which have not the same dimensions.
# My initial problem
xyplot(y ~ time, groups = char, data = foo,
panel = function(x, y, groups=groups, model = model, ...){
panel.xyplot(x,y, groups = groups, ...);
newdata <- data.frame(time = x, char = groups);
newdata$y <- predict(model, newdata = newdata);
# Here is an approximation of my results, which are good, except that
# I would like to increase the number of points for the curves
representation
# to obtain smooth representation.
panel.superpose(x = newdata$time, y=newdata$y, groups = groups, ...,
panel.groups = function(x,y, col, col.symbol, ...){
panel.lines(x,y, col.line = col.symbol)
})
}, model = res_nls);
#=============================
Many thanks,
Francois
2015-08-06 4:23 GMT+00:00 Duncan Mackay <dulcalma at bigpond.com>:
> Sending again as I pressed the wrong button
>
> Hi
>
> Your model is overfitted for the number of points - the fitted values bear
> no resemblance to a line fitting the data - ?too many x values to properly
> predict
> ? splines etc
> library(locfit)
> ?locfit
> If you want coefs of sort.
>
> You have to do a nls on each group to get the values.
> You have 2 ways to do this
> 1 do nls and predict y values separately as below an plot them in a panel
> function
> 2 do nls and predict within the panel function
>
> > res_nlsA <-
> + nls(y ~ SSgompertz(time, Asym = 100, b2, b3), data = subset(foo, char ==
> "A"),
> + start = list( b2 = c(7,7), b3 = c(0.9, 0.9)))
> There were 50 or more warnings (use warnings() to see the first 50)
> > res_nlsA
> Nonlinear regression model
> model: y ~ SSgompertz(time, Asym = 100, b2, b3)
> data: subset(foo, char == "A")
> b21 b22 b31 b32
> 4.6761 7.6280 0.9111 0.9000
> residual sum-of-squares: 10.59
>
> Number of iterations to convergence: 32
> Achieved convergence tolerance: 8.001e-06
> > res_nlsB <-
> + nls(y ~ SSgompertz(time, Asym = 100, b2, b3), data = subset(foo, char ==
> "B"),
> + start = list( b2 = c(7,7), b3 = c(0.9, 0.9)))
> > res_nlsB
> Nonlinear regression model
> model: y ~ SSgompertz(time, Asym = 100, b2, b3)
> data: subset(foo, char == "B")
> b21 b22 b31 b32
> 9.5231 86.6421 0.8164 0.6582
> residual sum-of-squares: 36.26
>
> Number of iterations to convergence: 10
> Achieved convergence tolerance: 1.37e-06
> >
> > predict(res_nlsB, newdata = list(time = seq(5,60,2.5)))
> [1] 3.162527 2.321040 28.575719 62.814774 63.489626 94.416582
> 84.809570 99.292612 94.199502 99.912322 97.856129 99.989162 99.217093
> 99.998661 99.715346
> [16] 99.999835 99.896669 99.999980 99.962512 99.999997 99.986402
> 100.000000 99.995068
> Warning messages:
> 1: In b3^x :
> longer object length is not a multiple of shorter object length
> 2: In -b2 * .expr2 :
> longer object length is not a multiple of shorter object length
> > predict(res_nlsB, newdata = list(time = seq(5,60,5)))
> [1] 3.162527 26.638924 63.489626 98.000694 94.199502 99.969171
> 99.217093 99.999529 99.896669 99.999993 99.986402 100.000000
> > predict(res_nlsA, newdata = list(time = seq(0,55,5)))
> [1] 0.9314904 1.1062577 15.8281748 20.7951276 48.3512870 57.8358643
> 75.0914782 82.6202586 89.3216331 93.5601696 95.6459638 97.7058236
>
> xyplot(y ~ time, data = foo,
> groups = char,
> panel = panel.superpose,
> panel.groups = function(x, y, ...,group.number){
>
> panel.xyplot(x,y, ...);
>
> xy.nls <-
> nls(y ~ SSgompertz(x, Asym = 100, b2, b3),
> start = list( b2 =7, b3 = 0.9))
>
> #print(summary(xy.nls))
>
> if (group.number == 1){
> xy.pre <- predict(xy.nls, newdata = list(time =
> seq(0,55,5)))
> print(xy.pre)
> llines( seq(5,60,5), xy.pre, col = "blue")
> } else {
> xy.pre <- predict(xy.nls, newdata = list(time =
> seq(5,60,5)))
> print(xy.pre)
> llines( seq(5,60,5), xy.pre, col = "magenta")
> }
>
> })
>
> Regards
>
> Duncan
>
> Duncan Mackay
> Department of Agronomy and Soil Science
> University of New England
> Armidale NSW 2351
> Email: home: mackay at northnet.com.au
>
>
> -----Original Message-----
> From: R-help [mailto:r-help-bounces at r-project.org] On Behalf Of François
> Collin
> Sent: Wednesday, 5 August 2015 22:24
> To: r-help at r-project.org
> Subject: [R] Lattice: raw data and prediction of a non linear fitted
> function
>
> Dear all,
>
> I have a question about lattice use. I would like to use it to represent:
> - my raw data as points,
> - and the results of a non linear fit as a line,
> - having 2 groups data (and so 2 colors).
>
> However, as I have few raw data, I would like to increase the number of
> points to smooth the line which correspond to the fit.
>
> So my questions are:
> - How can I use the group argument to make new predictions inside panel
> argument ?
> - How can I use this prediction inside the panel argument to draw the fit
> per group?
>
> Hereafter a minimal example:
>
> #==================================================
> library(lattice)
> set.seed(2)
>
> # Dummy dataframe
> foo <- data.frame(
> time = seq(0, 60, 5),
> char = sample(c("A", "B"), size = 13, replace = TRUE)
> );
>
> # Simulated response vector according a Gompertz function + rnorn(0, 5)
> foo$y <- SSgompertz(foo$time, Asym = 100, b2 = ifelse(foo$char == 'A', 6,
> 10),
> b3 = ifelse(foo$char == "A", .91, .8)) +
> rnorm(nrow(foo), mean=0, sd = 5);
>
> # Non-linear fit on simulation data
> res_nls <- nls(
> y ~ SSgompertz(time, Asym = 100, b2[char], b3[char]), data = foo,
> start = list( b2 = c(7,7), b3 = c(0.9, 0.9)));
>
> # My problem
> xyplot(y ~ time, groups = char, data = foo,
> panel = function(x, y, groups=groups, model = model, ...){
> panel.xyplot(x,y, groups = groups, ...);
>
> newdata <- data.frame(time = x, char = groups);
> newdata$y <- predict(model, newdata = newdata);
>
> panel.superpose(x = newdata$time, y=newdata$y, groups = groups, ...,
> panel.groups = function(x,y, col, col.symbol, ...){
> panel.lines(x,y, col.line = col.symbol)
> })
>
> }, model = res_nls);
> #==================================================
>
>
> Many thanks,
> Francois
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.
>
>
[[alternative HTML version deleted]]
More information about the R-help
mailing list