This lecture covers the following section in An Introduction to R:
One of the nice features of R is its graphics subroutines. Let us starts with the time series econ.dat
> X11() # To create a graphics window > x<-scan(file="econ.dat",skip=2) > ts.plot(x)To get hard copies, use (pp. 84-85)
> postscript("econ.ps", horizontal=F,height=6.width=4)
> ts.plot(x)
> dev.off()
The graph is in the file "econ.ps".
Histogram (related to density) and sample distribution function can be plotted by (pp. 41-42).
> par(mfrow=c(2,2)) > hist(x) # default output > hist(x,nclass=20) # with 20 interval > hist(x,nclass=20,probability=T) # with density label in the y-axis > library(stepfun) > plot(ecdf(x)) # CDF
Note that the par(mfrow=c(2,2)) instruction is the partition one page into a (2,2) compartments. It can be (m,n) for any m and n. Here are plots of some of the density functions (p.40, the preface d is for density)
postscript("density.ps")
par(mfrow=c(2,2))
x<-seq(-4,4,length=200)
y<-dnorm(x)
plot(x,y,type="l");title("Normal density")
x<-seq(0,10,length=200)
y<-dexp(x)
plot(x,y,type="l");title("Exponential density")
x<-seq(0,10,length=200)
y<-dgamma(x,2)
plot(x,y,type="l");title("Gamma(2) density")
x<-seq(0,10,length=200)
y<-dlnorm(x)
plot(x,y,type="l");title("Log-normal density")
dev.off()
Dr. Rajiv Patel's Data, T Raji Molt 40 0.840 1.150 38 0.930 1.150 36 1.050 1.210 34 1.020 1.130 32 1.050 1.180 30 1.150 1.290 28 1.140 1.350 26 1.250 1.430 24 1.290 1.530 22 1.410 1.660 20 1.530 1.720 18 1.590 1.770 16 1.680 1.930 14 1.820 2.070 12 1.890 2.180 10 2.020 2.240 8 2.020 2.380 6 2.220 2.510 4 2.290 2.670The following program is used to view the data.
postscript("patel.ps",height=5,width=6)
ex1<-scan("patel.dat",list(0,0,0),skip=2)
TC<-ex1[[1]]
RAJ<-ex1[[2]]
MOLT<-ex1[[3]]
x<-c(2:42) # producing theoretically fitted curves
cv1<- 3.0105*10**19/(273+x)**7.82529
cv2<- 3.46524*10**19/(273+x)**7.82529
x0<-5
y01<-0.8
y02<-0.5
x1<-c(8,14)
y1<-c(0.8,0.8)
y2<-c(0.5,0.5)
plot(TC,RAJ,xlim=range(0,50),ylim=range(0,3.5),
xlab="temperature (degree C)",
ylab="apparent microviscosity",
type="p",pch=0) # pch =0 means empty square
# Add more points (overlay)
points(TC,MOLT,type="p",pch=18) # pch =18 means black diamond
title("Raji and Molt microviscosity plot")
lines(x,cv1) # Add a straight line
lines(x,cv2,lty=4) # add another line
points(x0,y01,type="p",pch=0) # add a point at (x0,y01)
points(x0,y02,type="p",pch=18)
lines(x1,y1)
lines(x1,y2,lty=4)
text(19,0.8,"Raji") # text "Raji" at (19, 0.8)
text(19,0.5,"Molt")
dev.off()
If you wish to fit with a straight (linear regression) line. The computer
will find the regression line for you.
postscript("linearfit.ps",height=4,width=5)
plot(MOLT,RAJ,xlim=range(0,3),ylim=range(0,3),
xlab="MOLT viscosity", ylab="RAJ viscosity",type="p",pch=3) # pch 3 = +
abline(lsfit(MOLT,RAJ),lty=2) # lsfit=least square fit
# lty (line type) = = => dashed line
dev.off()
Longitudinal data display (high-level overlay)
The following data were the weights over time (weeks) of rats under the treatemnt of thiourocil versus control. Both groups had 10 individuals. It is saved as longit.dat.
1 2 3 4 5
CON 57 86 114 139 172
CON 60 93 123 146 177
CON 52 77 111 144 185
CON 49 67 100 129 164
..................
THI 56 78 95 103 108
THI 58 69 93 116 140
THI 46 61 78 90 107
THI 53 72 89 104 122
To disply them, we may use the following.
postscript("longit.ps",height=6,width=7.5)
x<-scan("longit.dat",list("",0,0,0,0,0),skip=1)
ym<-matrix(0,20,5)
for (i in c(1:5)) ym[,i]<-x[[i+1]]
ym # To make sure that this ym recovers the original matrix
xaxis<-c(1:5)
plot(xaxis,ym[1,],xlim=range(0,6),ylim=range(50,180), xlab="time in weeks",
ylab="Weight",type="both",lty=1,col="blue")
# type="both" means both points "p" and line "l"
for (i in c(2:10)) points(xaxis,ym[i, ],type="both",lty=1,col="blue")
for (i in c(11:20)) points(xaxis,ym[i, ],type="both",lty=1,col="red")
title("Thyroxin (red) vs control (blue)")
dev.off()
Multivariate data display
If we want to know the time relation between the weights, we can use
pairs() (p. 73). (Apparently, the relation losened when times were far away.)
postscript("pairs.ps")
pairs(ym)
dev.off()
Three dimensional graphics
Here is a 3-D graph using persp() (p.74).
postscript("twoballs.ps",height=6,width=6)
x<-seq(0,1,length=100)
y<-seq(0,1,length=100)
z<-matrix(0,100,100)
for(i in 1:100){for(j in 1:100){if ((x[i]-.5)^2+(y[j]-0.5)^2<0.09) z[i,j]<-(0.09-(x[i]-0.5)^2-(y[j]-0.5)^2)^0.5}}
for(i in 1:100){for(j in 1:100){if ((x[i]-.75)^2+(y[j]-0.25)^2<0.04) z[i,j]<-(0.04-(x[i]-0.75)^2-(y[j]-0.25)^2)^0.5}}
persp(x=seq(0,1,len=100),y=seq(0,1,len=100),z,xlim=range(0,1),
ylim=range(0,1),zlim=range(0,1.0))
dev.off()
