plot(resid(ldh.out1),type="l")
lines(resid(ldh.out1)[sort.list(ldh.df$v)],type="l",lty=2,col=2)
lines(resid(ldh.out1)[sort.list(ldh.df$a)],type="l",lty=4,col=4)
lines(resid(ldh.out1)[sort.list(ldh.df$b)],type="l",lty=6,col=6)
lines(resid(ldh.out1)[sort.list(ldh.df$p)],type="l",lty=8,col=8)
summary(ldh.out1)$parameters
summary(ldh.out1)$df
1-pt(2.554559,20)
1-pt(min(summary(ldh.out1)$parameters[,3]),summary(ldh.out1)$df[2])
?pt
min(summary(ldh.out1)$parameters[,3])
summary(ldh.out1)$df[2]
summary(ldh.out1)$parameters
1-pt(2.554559,20)
1-pnorm(2.554559)
eps<-qt((1 - 67/100)/2,20)*summary(ldh.out1)$parameters[1,2]
summary(ldh.out1)$parameters[1,1] + eps
summary(ldh.out1)$parameters[1,1] - eps
eps/summary(ldh.out1)$parameters[1,2]
qt((1 - 67/100)/2,20)
(1 - 67/100)/2
summary(ldh.out1)$parameters[1,2]
eps
ls()
rm(eps)
     qt(.975, df = c(1:10,20,50,100,1000))
     1 - pt(1:5, df = 1)
attach(ldh.df)
range(v0)
mean(v0)
summary(ldh.out1)$sigma
x<-c(min(v0), mean(v0), max(v0))
cbind(x,summary(ldh.out1)$sigma/x*100)
detach(2)
rm(x)
ls()
ls(pos=2)
attach(ldh.df)
range(v0)
mean(v0)
summary(ldh.out1)$sigma
x<-c(min(v0), mean(v0), max(v0))
cbind(x,summary(ldh.out1)$sigma/x*100)
detach(2)
rm(x)
q()
ls()
x11()
ldh.plot.df<-as.data.frame(cbind(ldh.df,v0calc=(ldh.df$v0 - ldh.out1$res)))
attach(ldh.plot.df)
ldh.df$v0
ldh.out1$res
ldh.out1$residuals
residuals(ldh.out1)
library()
?library
b
packages()
.packages()
(.packages())
plot.nls(ldh.out1)
plot.lme(ldh.out1)
?plot.lme
help.search("plot.lme"
)
library(nlme)
plot.lme(ldh.out1)
plot.nls(ldh.out1)
plot.lme(ldh.out1, v0 - fitted(.))
ldh.out1
attach(ldh.out1)
plot.lme(ldh.out1, v0 - fitted(.))
plot.lme(ldh.out1, .v0 - fitted(.))
ldh.out1$v0
ldh.out1.data
ldh.out1[data]
v0
as.list(ldh.out1)
list(ldh.out1)
?nlsModel
ldh
ls()
attach ldh.df
attach( ldh.df)
ls(0
)
ls()
ls(,2)
plot.lme(ldh.out1, v0 - fitted(.))
fitted(ldh.out1)
plot.lme(ldh.out1, v0 )
plot.lme(ldh.out1, v0 - residuals(.))
plot.lme(ldh.out1, resid(.) - fitted(.))
attach(ldh.out1)
ls(,2)
plot.lme(ldh.out1, resid(.) - fitted(.))
list(call)
list(data)
list(dataClasses)
list(m)
m$resid
m$resid()
ldh.out1$resid()
ldh.out1[m]$resid()
ldh.out1[[m]]$resid()
ldh.out1$$m$resid()
m$Rmat()
ldh.out1[m]$resid()
m$fitted()
ldh.df$v0
ls()
ls(,2)
 resdid(ldh.out1)
 resid(ldh.out1)
ldh.plot.df<-as.data.frame(cbind(ldh.df,v0calc=(ldh.df$v0 - resid(ldh.out1))))
ldh.plot.df
attach(ldh.plot.df)
#               For scatterplots
pairs( ~ v0 + a + b + p)
ldh.plot.df$p==0
pairs( ~ I(1/v0calc) + I(1/a) + I(1/b) + p)
nop<-ldh.plot.df[ldh.plot.df$p==0,]
nop
plot(1/nop$a,1/nop$v0)
points(1/nop$a,1/nop$v0calc,pch="+")
points(1/nop$b,1/nop$v0calc,pch="+")
points(1/nop$b,1/nop$v0calc,pch="+")
plot(1/nop$b,1/nop$v0)
points(1/nop$b,1/nop$v0calc,pch="+")
plot(ldh.plot.df)
?plot
plot.data.frame(ldh.plot.df)
methods(plot)
plot.histogram(ldh.plot.df)
plot.lme(ldh.plot.df)
plot.default(ldh.plot.df)
plot.ecdf(ldh.plot.df)
plot(nop)
plot(ppoints(v0), sort(v0))
plot(ppoints(v0-v0calc), sort(v0-v0calc))
ls(,2)
vo
v0
detach()
detach(2)
ls(,2)
ls(,2)
detach()
ls(,2)
detach()
ls(,2)
ls(,2)
attach(ldh.plot.df)
ls(,2)
plot(ppoints(v0-v0calc), sort(v0-v0calc))
plot(ppoints(v0), sort(v0))
plot(ppoints(v0-v0calc), sort(v0-v0calc))
qqnorm(v0-v0calc)
qqline(v0-v0calc)
boxplot(v0,v0calc)
boxplot(v0-v0calc)
nop<-ldh.plot.df[ldh.plot.df$p==0,]
nop.uniq <- nop[!duplicated(paste(nop[, "a"], nop[, "b"])),]
nop.uniq
nop
ls()
ls(,pos=2)
detach()
ls()
detach()
ls(,2)
?library
(.packages())
library("akima")
?interp
?persp
(.packages())
library("akima")
library("graphics")
(.packages())
i<-interp(1/nop.uniq$a, 1/nop.uniq$b, 1/nop.uniq$v0)
j<-interp(1/nop.uniq$a, 1/nop.uniq$b, 1/nop.uniq$v0calc)
k<-interp(1/nop.uniq$a, 1/nop.uniq$b, 1/nop.uniq$v0 - 1/nop.uniq$v0calc)
l<-interp(nop.uniq$a, nop.uniq$b, nop.uniq$v0)
m<-interp(nop.uniq$a, nop.uniq$b, nop.uniq$v0calc)
n<-interp(nop.uniq$a, nop.uniq$b, nop.uniq$v0 - nop.uniq$v0calc)
ls()
persp(i)
persp(j)
persp(k)
persp(l)
persp(m)
persp(n)
?contour
contour(i,nlevels=30)
contour(j,nlevels=30)
contour(k,nlevels=30)
contour(l,nlevels=30)
contour(m,nlevels=30)
contour(n,nlevels=30)
detach(2)
ls(,2)
ls()
rm(??????)
xx<-0
xx<ls()
rm(list=xx)
xx<-ls()
rm(list=xx)
ls)_
ls()
(.packages())
ls(,2)
ldh.df<-read.table("ldh.data",header=T)
#
# list data frame, extract a matrix element
ldh.df
ldh.df[1,3]             # gives an error?
as.matrix(ldh.df)[1,3]
#
# A data frame is _not_ a matrix.  (It is a list.)
as.list(ldh.df)
parameters.ldh<-list(Vmax=1e+00,KmA=1e-03,KmB=1e-02,KmAB=1e-03,C1=1e-03,C2=1e-06)
# 
# View the results.
parameters.ldh
# 
# Read the documentation on nls().
?nls
ldh.out1<-nls(v0 ~ Vmax/(1 + C1*p + KmA/a + (KmB/b)*(1 + C1*p)*(1 + C2*p) +
(KmAB/(a*b))*(1 + C1*p)), ldh.df, start=parameters.ldh, trace=TRUE)
#
#look at the results
as.list(ldh.out1)
summary(ldh.out1)
residuals(ldh.out1)
fitted.values(ldh.out1)
ldh.df$v0-fitted.values(ldh.out1) - residuals(ldh.out1)
coef(ldh.out1)
summary(ldh.out1)$parameters
summary(ldh.out1)$sigma
sqrt(sum(residuals(ldh.out1)^2)/summary(ldh.out1)$df[2])
ls()
ls(,2)
vcov(ldh.out1)
#(for S following works for var/covar matrix, but not needed with R)
#rinv<-backsolve(ldh.out1$R,diag(summary(ldh.out1)$df[1]));rinv %*% t(rinv)
#
# starting estimates need not be names - could be vector elements
nls(v0 ~ C[1]/(1 + C[5]*p + C[2]/a + (C[3]/b)*(1 + C[5]*p)*(1 + C[6]*p) + (C[4]/(a*b))*(1 + C[5]*p)), data=ldh.df,start=list(C=c(1e+00,1e-03,1e-02,1e-03,1e-03,1e-06)))
# cleanup
# rm(rinv) # not with R
ldh.out1
ls()
ttach(ldh.df)
objects(2)
ldh.nop.df<-data.frame(v0,a,b)[1:20,]
ls()
ls(,2)
ls()
attach(ldh.df)
ls()
ls(,2)
objects(2)
detach(2)
ls(,2)
attach(ldh.df)
objects(2)
ldh.nop.df<-data.frame(v0,a,b)[1:20,]
ldh.nop.df
ldh.df
ldh.nop.out1<-nls(v0 ~ (Vmax/(1 + KmA/a + KmB/b + KmAB/(a*b))), ldh.nop.df,list("Vmax"=1e+00,"KmA"=1e-03,"KmB"=1e-02,"KmAB"=1e-03))
#               And now a linear fit of the reciprocal rate equation
#               with lm().
lm.out1<-lm(1/v0 ~ I(1/a)*I(1/b),data=ldh.nop.df)
#               And weight the linear fit by the theoretical v0^4.
lm.wt.v04.out1<-lm(1/v0 ~ I(1/a)*I(1/b),data=ldh.nop.df, weights=ldh.nop.df$v0^4)
x<-coef(lm.wt.v04.out1)
(coef(ldh.nop.out1) - x/c(x[1]^2,x[1],x[1],x[1]))/coef(ldh.nop.out1)
x<-coef(lm.out1)
coef(ldh.nop.out1) - x/c(x[1]^2,x[1],x[1],x[1]))/coef(ldh.nop.out1)
(coef(ldh.nop.out1) - x/c(x[1]^2,x[1],x[1],x[1]))/coef(ldh.nop.out1))
x
coef(lm.wt.v04.out1
)
x[1]
coef(ldh.nop.out1)
(coef(ldh.nop.out1) - x/c(x[1]^2,x[1],x[1],x[1]))/coef(ldh.nop.out1)
x<-coef(lm.out1)
(coef(ldh.nop.out1) - x/c(x[1]^2,x[1],x[1],x[1]))/coef(ldh.nop.out1)
c(x[1]^2,x[1],x[1],x[1])
x
x/c(x[1]^2,x[1],x[1],x[1])
coef(ldh.nop.out1
)
summary(lm.wt.v04.out1)$cov
gmasq<-sum(summary(lmout)$residuals^2)/16 
gmasq<-sum(summary(lm.wt.v04.out1)$residuals^2)/16 
sigmasq<-sum(summary(lm.wt.v04.out1)$residuals^2)/16 
summary(lm.wt.v04.out1)$cov*sigmasq
?covariance
?cov
cov(lm.wt.v04.out1)
?nlsModel
vcov(lm.wt.v04.out1)
summary(lm.wt.v04.out1)$cov*sigmasq
lm.to.nls <- function(lmout) {
#for linear fit to LDH data
#compute values and std.error of Vmax and KmA, etc.
#from values and std.errors of 1/Vmax and KmA/Vmax, etc.
#
#value: data frame with estimate, std.error and t value for each parameter
#argument: data frame returned by lm()
#
        std.error<-numeric(4)
        parameters<-numeric(4)
        sigmasq<-sum(summary(lmout)$residuals^2)/16
        cij<-summary(lmout)$cov*sigmasq    # same as vcov(lmout)
        x<-coef(lmout)
        parameters[1]<-1/x[1]
        std.error[1]<-sqrt(1/x[1]^4*cij[1,1])
        for (i in 2:4) {
                parameters[i]<-x[i]/x[1]
                std.error[i]<-sqrt((1/x[1])^2*cij[i,i] + (-x[i]/(x[1]^2))^2*cij[1,1] + 2*(1/x[1])*(-x[i]/(x[1]^2))*cij[1,i])
        }
        data.frame(parameters = parameters, std.error = std.error, t = parameters/std.error)
}
lm.to.nls
ls()
rm(gmasq)
ls()
rm(sigmasq, x)
ls()
parameters.ldh
x<-lm.to.nls(lm.wt.v04.out1)
x
summary(ldh.nop.out1)$parameters
(summary(ldh.nop.out1)$parameters - x)/x
summary(aov(1/v0 ~ I(1/a) + I(1/b) + I(1/a):I(1/b),data=ldh.nop.df, weights=ldh.nop.df$v0^4))
lm.wt.nokmab.out1<-lm(1/v0 ~ I(1/a) + I(1/b),data=ldh.nop.df, weights=ldh.nop.df$v0^4)
ldh.nop.plot.df<-as.data.frame(cbind(ldh.nop.df,v0calc=(ldh.nop.df$v0 - ldh.nop.out1$res)))
ldh.nop.df$v0
ldh.nop.out1$res
ldh.nop.plot.df<-as.data.frame(cbind(ldh.nop.df,v0calc=(ldh.nop.df$v0 - residuals(ldh.nop.out))))
ldh.nop.plot.df<-as.data.frame(cbind(ldh.nop.df,v0calc=(ldh.nop.df$v0 - residuals(ldh.nop.out1))))
lm.nokmab.plot.df<-as.data.frame(cbind(ldh.nop.df, v0calc=fitted(lm.wt.nokmab.out1)))
q()
?help.search
?nls
(.packages())
?plot
ls(,2)
?as
 
?as
?as.matrix
(.packages())
q()
?dget
ls()
ldh.out1.sgi<-dget("ldh.out1.sgi")
?dump
?dget
dget
restore
dump
?restore
?dump
?save
history()
ldh.simpfit.out.sgi<-dget("ldh.simpfit.out.sgi")
ldh.simpfit.out.sgi
ldh.simpfit.out.sgi$coef
list(ldh.simpfit.out.sgi)
ldh.simpfit.out.sgi$coef
ldh.out1$coef
coef(ldh.out1)
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
?x11
 
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
?axes
help.search("axes")
?axis
q()
?axis
axes
?axes
title
?title
?par
?par
?x11
ls()
history()
source("ldhplot.all.s")
ls()
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
?axis
?par
?par
?plot
?title
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
graphics.off()
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef,"file")
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef,"file")
?postscript
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef,"file")
?postscript
 ?dev.off
?nop
help.search("nop")
 ?dev.off
dev.set()
graphics.off()
dev.set()
graphics.off()
dev.cur()(
)
dev.cur()
device()
?device
is.null(device)
is.null(dev.cur())
postscript()
dev.cur()
dec.off()
dev.off()
x11()
dev.cur()
dev.cur() == "x11"
dev.cur() == 2
?name
?names
history()
names(dev.cur())
names(dev.cur()) == "X11"
dev.cur() == "X11"
postscript()
dev.cur()(
)
dev.cur()
dev.list()
graphics.off()
dev.list()
dev.cur()
cat ""
cat("")
?cat
help.search("operation")
help.search("command")
help.search("commands")
help.search("continue")
help.search("ignore")
help.search("dummy")
help.search("syntax")
?Syntax
?Control
?Control
next
continue
break
.
""
"."
T
cat("")
?postscript
source("ldhplot.all.s")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef,"file")
ls()
dev.cur()
simpfit.out.sgi
ldh.simpfit.out.sgi<-dget("ldh.simpfit.out.sgi")
ls()
ldh.out1.sgi<-dget("ldh.out1.sgi")
SIMPFIT.temp.sgi<-dget("SIMPFIT.temp.sgi")
dget
?dget
?structure
ls()
SIMPFIT.temp.sgi
l
?dget
dget("ldh.out1.sgi")
?match.call
dget
q()
ls()
ldh.simpfit.out.sgi$coef
l
ldh.simpfit.out.sgi$model
?system
 ?!
help("!")
 q()
source("Plots/ldhplot.all.s")
ldh.simpfit.out.sgi<-dget("Plots/ldh.simpfit.out.sgi")
ldhplot(ldh.simpfit.out.sgi$data, ldh.simpfit.out.sgi$coef)
history()
graphics.off()
q()
