[R] Plotting observed vs. fitted values

[Gray Calhoun]

oscar linares wrote:
> Dear Wiza[R]ds,
> 
> I am very grateful to Duncan Murdoch for his assistance with this problem.
> His help was invaluable. However, the problem has become a little more
> complicated for me. Now, in each plot, I need to plot the observed and
> fitted values of a supine and upright posture experiment. Here is what I
> have and how far I got.
> 
> # tritiated (3H)-Norepinephrine(NE) disappearance from plasma
> # concentrations supine and upright
> # supine
> datasu <- data.frame(
> time=c(0,1,2,4,6,8,10,15,20),
> concsu=c(385.61,265.88,173.87,99.47,66.7,55.27,48.29,39.85,40.66)
> )
> # upright
> dataup <- data.frame(
> time=c(0,1,2,4,6,8,10,15,20),
> concup=c(767.27,529.03,328.13,225.94,164,151.1,132.02,121.15,70.88)
> )
> 
> # fit supine
> wt.su<- function(resp, time,A1,a1,A2,a2)
> {
>     pred <- A1*exp(-a1*time)+A2*exp(-a2*time)
>     (resp - pred) / sqrt(pred)
> }
> 
> fit.wt.su <- nls( ~ wt.su(concsu, time,A1,a1,A2,a2), data=datasu,
>               start=list(A1=500,a1=1,A2=100,a2=0.2),
>               trace = TRUE)
> 
> summary(fit.wt.su)
> 
> # fit upright
> wt.up<- function(resp, time,B1,b1,B2,b2)
> {
>     pred <- B1*exp(-b1*time)+B2*exp(-b2*time)
>     (resp - pred) / sqrt(pred)
> }
> 
> fit.wt.up <- nls( ~ wt.up(concup, time,B1,b1,B2,b2), data=dataup,
>               start=list(B1=500,b1=1,B2=100,b2=0.1),
>               trace = TRUE)
> 
> summary(fit.wt.up)
> 
> # Function that returns predicted values and graphics
> # by Duncan Murdoch
> 
> predictionssu <- function(fit, time) {
>   params <- summary(fit)$coefficients[, 1]
>   A1 <- params["A1"]
>   a1 <- params["a1"]
>   A2 <- params["A2"]
>   a2 <- params["a2"]
> 
>   A1*exp(-a1*time)+A2*exp(-a2*time)
> }
> 
> predictionsup <- function(fit, time) {
>   params <- summary(fit)$coefficients[, 1]
>   B1 <- params["B1"]
>   b1 <- params["b1"]
>   B2 <- params["B2"]
>   b2 <- params["b2"]
> 
>   B1*exp(-b1*time)+B2*exp(-b2*time)
> }
> 
> # plot observed and predicted values supine and upright in
> # each plot type (linear and smi-log)
> 
> # this does supine
> times <- seq(0,20, len=100)
> par(mfrow=c(2,1))
> plot(concsu ~ time, data=datasu)
> lines(times, predictionssu(fit.wt.su, times))
> plot(concsu ~ time, data=datasu, log="y")
> lines(times, predictionssu(fit.wt.su, times))
> 
> # Need upright in same plots (e.g., dashed line)
> 
> 
> Many thanks in advance...

Hi Oscar,
   The "lines" and "points" commands add to the current plot, so adding 
more points and dashed lines is straightforward (see the help page for 
"par" for more details on how to customize the line type).  The 
following code plots the supine and upright graphs together. (obviously 
it needs to go at the end of your code).

 > par(mfrow=c(2,1))
 > plot(c(datasu$concsu, dataup$concup) ~ c(datasu$time, dataup$time), 
xlab = "time", ylab = "conc")
 > lines(times, predictionssu(fit.wt.su, times))
 > lines(times, predictionsup(fit.wt.up, times), lty = "dashed")
 > plot(c(datasu$concsu, dataup$concup) ~ c(datasu$time, dataup$time), 
xlab = "time", ylab = "conc", log = "y")
 > lines(times, predictionssu(fit.wt.su, times))
 > lines(times, predictionsup(fit.wt.up, times), lty = "dashed")

The 'plot' command is ugly and could be improved; also both "conc" and 
"time" should be given better labels.  Since you have different data 
frames for "supine" and "upright", you don't need to distinguish the row 
names; you could have "conc" be the row in both data frames.  Even 
better, you could have "supine" vs. "upright" as a factor, and store all 
of your data in a single data frame; doing so would make it easier to 
use the lattice graphics package for this sort of problem.

-- 
Gray Calhoun

Assistant Professor of Economics, Iowa State University
http://www.econ.iastate.edu/~gcalhoun/

467 Heady Hall
Ames, IA 50011
Phone: (515) 294-6271
Fax:   (515) 294-0221