R_in_action

1672 days ago by takepwave

# ChromeではPDFをimgタグで表示できないので、pngに変換するように変更 2012/06/27 # # Rのグラフをsageで表示するためのユーティリティ関数 import time def preGraph(pdfFile): filename = DATA+pdfFile r.pdf(file='"%s"' %filename) return pdfFile def offGraph(): r.dev_off() def postGraph(pdfFile, fac=0.75): r.dev_off() width = int(640*fac) # html('<img src="%s?%s" width="%spx">'%(pdfFile, time.time(), width)) pngFile = convertPdf2Png(pdfFile) html('<img src="%s?%s" width="%spx">'%(pngFile, time.time(), width)) def getGraph(pdfFile, fac=0.75): width = int(640*fac) # return '<img src="%s?%s" width="%spx">'%(pdfFile, time.time(), width) pngFile = convertPdf2Png(pdfFile) return '<img src="%s?%s" width="%spx">'%(pngFile, time.time(), width) def convertPdf2Png(pdfFile): pngFile = pdfFile.replace(".pdf", ".png") result = os.popen("cd %s; convert -density 600 -geometry 1000 %s %s"%(DATA, pdfFile, pngFile)).read() # print result return pngFile 
       
# 変数の受け渡し # Sage -> R # すべてリストの形式に変換されてRに渡される s_age = r([1,3,5,2,11,9,3,9,12,3]).name('age') # name関数で名前を指定することで、Rの変数も定義される print r('age') # RからSageへの変換は、sageobj関数または_sage_()メソッドを使用 print sageobj(r('age')) 
       
 [1]  1  3  5  2 11  9  3  9 12  3
[1, 3, 5, 2, 11, 9, 3, 9, 12, 3]
 [1]  1  3  5  2 11  9  3  9 12  3
[1, 3, 5, 2, 11, 9, 3, 9, 12, 3]
# 文字列配列の受け渡しができないことが分かりました。以下のように無理矢理セットしてください。 ary = ["a", "b"] ss = str(ary)[1:-1]; print ss # リストを文字列に変換して[と]を取り除く rary = r('c(%s)'%ss).name('rary') print rary # Rに変換した配列 sary = sageobj(rary) print sary # sageに戻した配列 
       
'a', 'b'
[1] "a" "b"
['a', 'b']
'a', 'b'
[1] "a" "b"
['a', 'b']
# マトリックスの生成も同様に生成後に名前をつけて代入すればよい m = r.matrix([1, 2, 3, 4], 2, 2, byrow=True).name('m'); m 
       
     [,1] [,2]
[1,]    1    2
[2,]    3    4
     [,1] [,2]
[1,]    1    2
[2,]    3    4
# rでも同じ変数名が定義される r('m') 
       
     [,1] [,2]
[1,]    1    2
[2,]    3    4
     [,1] [,2]
[1,]    1    2
[2,]    3    4
# Rのマトリックス変数をsageの変数に変換する print type(m) sm = sageobj(m) show(sm) print type(sm) 
       
<class 'sage.interfaces.r.RElement'>

<type 'sage.matrix.matrix_integer_dense.Matrix_integer_dense'>
<class 'sage.interfaces.r.RElement'>

<type 'sage.matrix.matrix_integer_dense.Matrix_integer_dense'>
# fig1.4 age = r([1,3,5,2,11,9,3,9,12,3]).name('age') weight = r([4.4,5.3,7.2,5.2,8.5,7.3,6.0,10.4,10.2,6.1]).name('weight') print r.mean(weight) print r.sd(weight) print r.cor(age, weight) 
       
[1] 7.06
[1] 2.077498
[1] 0.9075655
[1] 7.06
[1] 2.077498
[1] 0.9075655
# グラフ表示 graph = preGraph("fig1.4.pdf") r.plot(age, weight) postGraph(graph, fac=0.5) 
       
# タイトルを追加する場合には、r.plotではなく、r関数を使用する graph = preGraph("fig1.4withTitle.pdf") r("plot(age, weight, main='Title')") postGraph(graph, fac=0.5) 
       
# fig 3.1 fig_3_1 = preGraph("fig3.1.pdf") r('attach(mtcars)') r('plot(wt, mpg)') r('abline(lm(mpg~wt))') r('title("Regression of MPG on Weight")') r('detach(mtcars)') postGraph(fig_3_1, fac=0.5) 
       
# fig 3.2 dose = r.c([20, 30, 40, 45, 60]).name('dose') drugA = r.c([16, 20, 27, 40, 60]).name('drugA') drugB = r.c([15, 18, 25, 31, 40]).name('drugB') 
       
fig_3_2 = preGraph("fig3.2.pdf") r('plot(dose, drugA, type="b")') postGraph(fig_3_2, fac=0.5) 
       
# fig 3.3 fig_3_3 = preGraph("fig3.3.pdf") opar = r('par(no.readlonly=TRUE)').name('opar') r('par(lty=2, pch=17)') r('plot(dose, drugA, type="b")') r.par(opar) postGraph(fig_3_3, fac=0.5) 
       
# fig 3.7 fig_3_7a = preGraph("fig3.7a.pdf") r('par(pin=c(2, 3))') r('par(lwd=2, cex=1.5)') r('par(cex.axis=.75, font.axis=3)') r('plot(dose, drugA, type="b", pch=19, lty=2, col="red")') offGraph() 
       
# parは、プロット毎にセットしないとダメ fig_3_7b = preGraph("fig3.7b.pdf") r('par(pin=c(2, 3))') r('par(lwd=2, cex=1.5)') r('par(cex.axis=.75, font.axis=3)') r('plot(dose, drugB, type="b", pch=23, lty=6, col="blue", bg="green")') r.par(opar) offGraph() html.table([[getGraph(fig_3_7a, fac=0.5), getGraph(fig_3_7b, fac=0.5)]]) 
       
# fig 3.8 fig_3_8 = preGraph("fig3.8.pdf") r('plot(dose, drugA, type="b", col="red", lty=2, pch=2, lwd=2, main="Clinical Trials for Drug A", sub="This is hypothetical data", xlab="Dosage", ylab="Drug Response", xlim=c(0, 60), ylim=c(0, 70))') postGraph(fig_3_8, fac=0.5) 
       
# fig 3.9 x = r('c(1:10)').name('x') y = r('x').name('y') z = r('10/x').name('z') print x print y print z 
       
 [1]  1  2  3  4  5  6  7  8  9 10
 [1]  1  2  3  4  5  6  7  8  9 10
 [1] 10.000000  5.000000  3.333333  2.500000  2.000000  1.666667 
1.428571  1.250000  1.111111
[10]  1.000000
 [1]  1  2  3  4  5  6  7  8  9 10
 [1]  1  2  3  4  5  6  7  8  9 10
 [1] 10.000000  5.000000  3.333333  2.500000  2.000000  1.666667  1.428571  1.250000  1.111111
[10]  1.000000
fig_3_9 = preGraph("fig3.9.pdf") opar = r('par(no.readlonly=TRUE)').name('opar') r('par(mar=c(5, 4, 4, 8) + 0.1)') r('plot(x, y, type="b", pch=21, col="red", yaxt="n", lty=3, ann=FALSE)') r('lines(x, z, type="b", pch=22, col="blue", lty=2)') r('axis(2, at=x, labels=x, col.axis="red", las=2)') r('axis(4, at=z, labels=round(z, digits=2), col.axis="blue", las=2, cex.axis=0.7, tck=-.01)') r('mtext("y=1/x", side=4, line=3, cex.lab=1, las=2, col="blue")') r('title("An Example of Creative Axes", xlab="X values", ylab="Y=X")') r.par(opar) postGraph(fig_3_9, fac=0.5) 
       
# fig3.10 r('library(Hmisc)') fig_3_10 = preGraph("fig3.10.pdf") opar = r('par(no.readlonly=TRUE)').name('opar') r('par(lwd=2, cex=1.5, font.lab=2)') r('plot(dose, drugA, type="b", pch=15, lty=1, col="red", ylim=c(0, 60), main="Drug A vs. Drug B", xlab="Drug Dosage", ylab="Drug Response")') r('lines(dose, drugB, type="b", pch=17, lty=2, col="blue")') r('abline(h=c(30), lwd=1.5, lty=2, col="gray")') r('minor.tick(nx=3, ny=3, tick.ratio=0.5)') r('legend("topleft", inset=.05, title="Drug Type", c("A","B"), lty=c(1, 2), pch=c(15, 17), col=c("red", "blue"))') r.par(opar) postGraph(fig_3_10, fac=0.5) 
       
# fig 3.11 r('attach(mtcars)') fig_3_11 = preGraph("fig3.11.pdf") r('plot(wt, mpg, main="Mileage vs. Car Weight", xlab="Weight", ylab="Mileage", pch=18, col="blue")') r('text(wt, mpg, row.names(mtcars), cex=0.6, pos=4, col="red")') postGraph(fig_3_11, fac=0.7) r('detach(mtcars)') 
       

<environment: 0x30129c0>
attr(,"name")
[1] "mtcars"

<environment: 0x30129c0>
attr(,"name")
[1] "mtcars"
# fig 3.12 fig_3_12 = preGraph("fig3.12.pdf") opar = r('par(no.readlonly=TRUE)').name('opar') r('par(cex=1.5)') r('plot(1:7,1:7,type="n")') r('text(3,3,"Example of default text")') r('text(4,4,family="mono","Example of mono-spaced text")') r('text(5,5,family="serif","Example of serif text")') r.par(opar) postGraph(fig_3_12, fac=0.5) 
       
# fig 3.14 fig_3_14 = preGraph("fig3.14.pdf") r('attach(mtcars)') opar = r('par(no.readlonly=TRUE)').name('opar') r('par(mfrow=c(2,2))') r('plot(wt,mpg, main="Scatterplot of wt vs. mpg")') r('plot(wt,disp, main="Scatterplot of wt vs disp")') r('hist(wt, main="Histogram of wt")') r('boxplot(wt, main="Boxplot of wt")') r.par(opar) r('detach(mtcars)') postGraph(fig_3_14, fac=0.6) 
       
# fig 3.15 fig_3_15 = preGraph("fig3.15.pdf") r('attach(mtcars)') opar = r('par(no.readlonly=TRUE)').name('opar') r('par(mfrow=c(3,1))') r('hist(wt)') r('hist(mpg)') r('hist(disp)') r.par(opar) r('detach(mtcars)') postGraph(fig_3_15, fac=0.6) 
       
# fig 3.18 fig_3_18 = preGraph("fig3.18.pdf") opar = r('par(no.readlonly=TRUE)').name('opar') # Set up scatter plot r('par(fig=c(0, 0.8, 0, 0.8))') r('plot(mtcars$wt, mtcars$mpg, xlab="Miles Per Gallon", ylab="Car Weight")') # Add box plot above r('par(fig=c(0, 0.8, 0.55, 1), new=TRUE)') r('boxplot(mtcars$wt, horizontal=TRUE, axes=FALSE)') # Add box plot to right r('par(fig=c(0.65, 1, 0, 0.8), new=TRUE)') r('boxplot(mtcars$mpg, axes=FALSE)') r('mtext("Enhanced Scatterplot", side=3, outer=TRUE, line=-3)') r.par(opar) postGraph(fig_3_18, fac=0.6) 
       
# Sec 6.1.1 Simple bar plots # vcdパッケージを使用するため、1回のみ実行 # r('install.packages("vcd")') 
       
r("library(vcd)") 
       
[1] "vcd"       "grid"      "stats"     "graphics"  "grDevices" "utils" 
"datasets"  "methods"  
[9] "base"     
[1] "vcd"       "grid"      "stats"     "graphics"  "grDevices" "utils"     "datasets"  "methods"  
[9] "base"     
counts = r("table(Arthritis$Improved)").name('counts'); counts 
       
  None   Some Marked 
    42     14     28 
  None   Some Marked 
    42     14     28 
# fig 6.1a fig6_1a = preGraph("fig6.1a.pdf") r('barplot(counts, main="Simple Bar Plot", xlab="Improvement", ylab="Frequency")') offGraph() 
       
# fig 6.1b fig6_1b = preGraph("fig6.1b.pdf") r('barplot(counts, main="Horizontal Bar Plot", xlab="Frequency", ylab="Improvement", horiz=TRUE)') offGraph() 
       
html.table([[getGraph(fig6_1a, fac=0.5), getGraph(fig6_1b, fac=0.5)]]) 
       
# Sec. 6.2 counts = r("table(Arthritis$Improved, Arthritis$Treatment)").name('counts'); counts 
       
        
         Placebo Treated
  None        29      13
  Some         7       7
  Marked       7      21
        
         Placebo Treated
  None        29      13
  Some         7       7
  Marked       7      21
# fig 6.2 fig6_2a = preGraph("fig6.2a.pdf") r('barplot(counts, main="Stacked Bar Plot", xlab="Treatment", ylab="Frequency", col=c("red", "yellow","green"), legend=rownames(counts))') offGraph() fig6_2b = preGraph("fig6.2b.pdf") r('barplot(counts, main="Grouped Bar Plot", xlab="Treatment", ylab="Frequency", col=c("red", "yellow", "green"), legend=rownames(counts), beside=TRUE)') offGraph() html.table([[getGraph(fig6_2a, fac=0.5), getGraph(fig6_2b, fac=0.5)]]) 
       
# fig 6.5 fig6_5 = preGraph("fig6.5.pdf") r('attach(Arthritis)') r('counts <- table(Treatment, Improved)') r('spine(counts, main="Spinogram Example")') r('detach(Arthritis)') postGraph(fig6_5, fac=0.6) 
       
# plotrixパッケージを使用するため、1回のみ実行 # r('install.packages("plotrix")') 
       
# fig 6.6 fig6_6 = preGraph("fig6.6.pdf") r('par(mfrow=c(2, 2))') r('slices <- c(10, 12,4, 16, 8)') r('lbls <- c("US", "UK", "Australia", "Germany", "France")') r('pie(slices, labels = lbls, main="Simple Pie Chart")') r('pct <- round(slices/sum(slices)*100)') r('lbls2 <- paste(lbls, " ", pct, "%", sep="")') r('pie(slices, labels=lbls2, col=rainbow(length(lbls2)), main="Pie Chart with Percentages")') r('library(plotrix)') r('pie3D(slices, labels=lbls,explode=0.1, main="3D Pie Chart ")') r('mytable <- table(state.region)') r('lbls3 <- paste(names(mytable), "\n", mytable, sep="")') r('pie(mytable, labels = lbls3, main="Pie Chart from a Table\n (with sample sizes)")') postGraph(fig6_6, fac=0.6) 
       
# fig 6.7 fig6_7 = preGraph("fig6.7.pdf") r('fan.plot(slices, labels = lbls, main="Fan Plot")') postGraph(fig6_7, fac=0.6) 
       
# fig 6.8 fig6_8 = preGraph("fig6.8.pdf") r('par(mfrow=c(2,2))') r('hist(mtcars$mpg)') r('hist(mtcars$mpg, breaks=12, col="red", xlab="Miles Per Gallon", main="Colored histogram with 12 bins")') r('hist(mtcars$mpg, freq=FALSE, breaks=12, col="red", xlab="Miles Per Gallon", main="Histogram, rug plot, density curve")') r('rug(jitter(mtcars$mpg))') r('lines(density(mtcars$mpg), col="blue", lwd=2)') r('x <- mtcars$mpg') r('h<-hist(x, breaks=12, col="red", xlab="Miles Per Gallon", main="Histogram with normal curve and box")') r('xfit<-seq(min(x), max(x), length=40)') r('yfit<-dnorm(xfit, mean=mean(x), sd=sd(x))') r('yfit <- yfit*diff(h$mids[1:2])*length(x)') r('lines(xfit, yfit, col="blue", lwd=2)') r('box()') postGraph(fig6_8, fac=0.7) 
       
# fig 6.9 fig6_9 = preGraph("fig6.9.pdf") r('par(mfrow=c(2,1))') r('d <- density(mtcars$mpg)') r('plot(d)') r('plot(d, main="Kernel Density of Miles Per Gallon")') r('polygon(d, col="red", border="blue")') r('rug(mtcars$mpg, col="brown")') postGraph(fig6_9, fac=0.6) 
       
# smパッケージを使用するため、1回のみ実行 # r('install.packages("sm")') 
       
# fig 6.10 # 凡例が表示されない fig6_10 = preGraph("fig6.10.pdf") r('par(lwd=2)') r('library(sm)') r('attach(mtcars)') r('cyl.f <- factor(cyl, levels= c(4,6,8), labels = c("4 cylinder", "6 cylinder", "8 cylinder"))') r('sm.density.compare(mpg, cyl, xlab="Miles Per Gallon")') r('title(main="MPG Distribution by Car Cylinders")') r('colfill<-c(2:(1+length(levels(cyl.f))))') r('legend(locator(1), levels(cyl.f), fill=colfill)') r('detach(mtcars)') postGraph(fig6_10, fac=0.6) 
       
# fig 6.14 fig6_14 = preGraph("fig6.14.pdf") r('mtcars$cyl.f <- factor(mtcars$cyl, levels=c(4,6,8), labels=c("4","6","8"))') r('mtcars$am.f <- factor(mtcars$am, levels=c(0,1), labels=c("auto", "standard"))') r('boxplot(mpg ~ am.f *cyl.f, data=mtcars, varwidth=TRUE, col=c("gold","darkgreen"), main="MPG Distribution by Auto Type", xlab="Auto Type")') postGraph(fig6_14, fac=0.6) 
       
# fig 6.16 fig6_16 = preGraph("fig6.16.pdf") r('dotchart(mtcars$mpg, labels=row.names(mtcars), cex=.7, main="Gas Mileage for Car Models", xlab="Miles Per Gallon")') postGraph(fig6_16, fac=0.6) 
       
# fig 11.1 fig11_1 = preGraph("fig11.1.pdf") r('attach(mtcars)') r('plot(wt, mpg, main="Basic Scatter plot of MPG vs. Weight", xlab="Car Weight (lbs/1000)", ylab="Miles Per Gallon ", pch=19)') r('abline(lm(mpg~wt), col="red", lwd=2, lty=1)') r('lines(lowess(wt,mpg), col="blue", lwd=2, lty=2)') postGraph(fig11_1, fac=0.6) 
       
# fig 11.3 fig11_3 = preGraph("fig11.3.pdf") r('pairs(~mpg+disp+drat+wt, data=mtcars, main="Basic Scatter Plot Matrix")') postGraph(fig11_3, fac=0.6) 
       
# 散布図がつぶれてしまったような場合の対処方法 r('set.seed(1234)') r('n <- 10000') r('c1 <- matrix(rnorm(n, mean=0, sd=.5), ncol=2)') r('c2 <- matrix(rnorm(n, mean=3, sd=2), ncol=2)') r('mydata <- rbind(c1, c2)') r('mydata <- as.data.frame(mydata)') r('names(mydata) <- c("x", "y")') 
       
[1] "x" "y"
[1] "x" "y"
fig11_7 = preGraph("fig11.7.pdf") r('with(mydata, plot(x, y, pch=19, main="Scatter Plot with 10,000 Observations"))') postGraph(fig11_7, fac=0.6) 
       
# fig 11.8 fig11_8 = preGraph("fig11.8.pdf") r('with(mydata, smoothScatter(x, y, main="Scatterplot Colored by Smoothed Densities"))') postGraph(fig11_8, fac=0.6) 
       
r('ftable(Titanic)') 
       
                   Survived  No Yes
Class Sex    Age                   
1st   Male   Child            0   5
             Adult          118  57
      Female Child            0   1
             Adult            4 140
2nd   Male   Child            0  11
             Adult          154  14
      Female Child            0  13
             Adult           13  80
3rd   Male   Child           35  13
             Adult          387  75
      Female Child           17  14
             Adult           89  76
Crew  Male   Child            0   0
             Adult          670 192
      Female Child            0   0
             Adult            3  20
                   Survived  No Yes
Class Sex    Age                   
1st   Male   Child            0   5
             Adult          118  57
      Female Child            0   1
             Adult            4 140
2nd   Male   Child            0  11
             Adult          154  14
      Female Child            0  13
             Adult           13  80
3rd   Male   Child           35  13
             Adult          387  75
      Female Child           17  14
             Adult           89  76
Crew  Male   Child            0   0
             Adult          670 192
      Female Child            0   0
             Adult            3  20
# fig 11.23 r('library(vcd)') fig11_23 = preGraph("fig11.23.pdf") r('mosaic(Titanic, shade=TRUE, legend=TRUE)') postGraph(fig11_23, fac=0.6)