Rはフリーの統計解析ソフトで、Sageのパッケージに付属しており、SageからRのコマンドを実行するための インターフェースも用意されています。
Rの特徴に以下のものがあげられます。
ここでは、SageからRのコマンドを使用する場合の注意点とRインタフェースの使用方法に説明します。
%%HTML
<link rel="stylesheet" type="text/css" href="css/sage_table_form.css">
SageからRのコマンドを使いやすくするために、以下のユーティリティ関数を用意、RUtil.pyに納めました。
load関数を使ってRUtil.pyを読み込みます。
# RとPandasのデータフレームを相互に変換する関数を読み込む
# Rの必要なライブラリ
r('library(ggplot2)')
r('library(jsonlite)')
# jupyter用のdisplayメソッド
from IPython.display import display, Latex, HTML, Math, JSON
# sageユーティリティ
load('script/sage_util.py')
# Rユーティリティ
load('script/RUtil.py')
Rとのインタフェースとして、r関数が提供されています。 r関数では、Rのコマンドを文字列で与えます。Rのコマンド内でダブルクォートを使っている場合には、シングルクォートでくくり、 Rのコマンド内でシングルクォートを使っている場合には、ダブルクォートを使って指定します。
r('Rのコマンド')RのコマンドにSageの変数の文字列を代入する場合には、以下のように%を使ってセットします。
r('Rのコマンド1%s, Rのコマンド2%s' % (変数1, 変数2))上記の例では、1番目の%sが変数1に、2番目の%sが変数2の値に置き換わります。
r関数の例として、Rで計算させた2のsqrtの結果を出力します。Rの出力では最初に[1]のように番号が付いて 表示されます。
# Rインターフェース
r('sqrt(2)')
[1] 1.414214
関数定義やRの制御コマンド(ループやif文等)を含む複数行のRコマンドをSageから実行する場合には、 ファイルにまとめてワークシートにアップロードして使います。
例でとして%%writefileコマンドを使って、data/mystats.txtに以下のように自乗を計算するx2関数を書いておきます。
%%writefile data/mystats.txt
x2 <- function(x) {
return (x^x)
}
Writing data/mystats.txt
ファイルに定義されたRのコマンドを読み込み評価するには、Rのsourceコマンドを使用します。
以下の例では、mystats.txtに書かれたx2の関数定義を実行し、リストlの各要素に関数x2を適応しています。
# Source関数でテキストファイルの内容を評価実行します
r('source("data/mystats.txt")')
# リストlの各要素に関数x2を適応する
r('l <- c(1, 2, 3)')
r('lapply(l, x2)')
[[1]] [1] 1 [[2]] [1] 4 [[3]] [1] 27
SageとRとの連携で重要なポイントはSageの変数をどうやってRに渡すかです。
以下の例では、変数slに[1, 2, 3]のリストをセットし、これをr関数に渡しています。 このままですとrlに返されるRElementオブジェクトの名前がSagexxのように自動で 割り当てられてしまします。
そこで、name関数でRにおける変数名を指定します。
rl = r(sl).name('l')説明でRの変数とSageの変数を混同しないようにRでの変数名をlとし、r関数の返した 値をrl変数にセットします。
出力結果から、r('l')で正しくslのリストがRに渡っていることが確認できます。 また、rlのタイプは、sage.interfaces.r.RElementであり、その値の出力結果は r('l')と同じであることも確認できます。
# Sageの変数(リスト)をRに渡す
sl = [1, 2, 3]
# RにSageの変数slを渡し、Rでの変数名をlとする
rl = r(sl).name('l')
print r('l')
print type(rl)
print rl
[1] 1 2 3 <class 'sage.interfaces.r.RElement'> [1] 1 2 3
リストの次にRでよく使われるのがマトリックスです。 例として、r.matrix関数を使って[1, 2, 3, 4]のリストから、2x2のマトリックスを作成します。 byrow=Trueとすることで、Sageと同じ行単位での配置を指定しています。 (RではFORTRANの影響で列単位で配列を保持しています)
# マトリックスの生成も同様に生成後に名前をつけて代入すればよい
m = r.matrix([1, 2, 3, 4], 2, 2, byrow=True).name('m'); m
[,1] [,2] [1,] 1 2 [2,] 3 4
Rに定義された変数mの値を出力すると、上記の出力と同じことであることが確認できます。
# rでも同じ変数名が定義される
r('m')
[,1] [,2] [1,] 1 2 [2,] 3 4
次にRのマトリックスをSageのマトリックスに変換する方法を示します。
先に定義されたRのマトリックスmをsageobj関数を使ってSageのマトリックスに変換します。 以下の例では、r('m')で返されたRElementをsageobjでSageの変数smに代入し、 その値とタイプを確認しています。smは正しく変換され、smのタイプも sage.matrix.matrix_integer_dense.Matrix_integer_denseとなっていることが 確認できます。
# 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'>
Sageの変数smからRのマトリックスに変換する場合には、listメソッドを使ってマトリックスsm をリストに変換し、smの行、sm列を使ってRのマトリックスに変換します。
# Sage => Rの場合には、listメソッドでマトリックスからリストに変換する
r.matrix(sm.list(), sm.nrows(), sm.ncols(), byrow=True)
[,1] [,2] [1,] 1 2 [2,] 3 4
残念ながら、Sageの文字列リストは、そのままではRの文字列リストに変換されません。 そこで、文字列リストをstr関数を使って文字列に変換し、最初[と最後の]を取り除いて、 r('c(%s)'%文字列リスト)を使ってRの文字列リストに変換します。
# 文字列のリストの場合には、上記の方法では上手くいかない
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']
Rのプロット関数を使う場合には、preGraphとpostGraphで挟むようにします。
# 正規分布をRインターフェースを使って表示
graph = preGraph("images/ch12_fig1.pdf")
r('x <- pretty(c(-3,3), 30)')
r('y <- dnorm(x)')
r('plot(x, y, type = "l", xlab = "Normal Deviate", ylab = "Density", yaxs = "i")')
postGraph(graph)
上記と同じことをRの関数dnormの値を使ってSageのplot関数でプロットした例を以下に示します。
# 同じグラフをRと連携して表示した例(計算にRの関数を使用)
# Sage 6.7ではN()を使った数値化が失敗するため、sageobjを使ってrの戻り値を数値に変換
x = var('x')
f = lambda x: sageobj(r.dnorm(x)) # r.dnormを使ってRのdnorm関数を呼び出す
plot(f, [x, -3, 3], axes_labels=["Normal Deviate", "Density"], figsize=6)
次に、r.plotを使ったプロットの例を示します。X軸にage、Y軸にweightを表示します。 r.plotにX, Yにsageの変数age, weightを使った場合、軸ラベル名にsagexxのような 名前が表示されます。
# グラフ表示(軸ラベルが正しく表示されていない)
age = [1,3,5,2,11,9,3,9,12,3]
weight = [4.4,5.3,7.2,5.2,8.5,7.3,6.0,10.4,10.2,6.1]
graph = preGraph("images/ch12_fig2.pdf")
r.plot(age, weight)
postGraph(graph)
age, weightにnameコマンドでRの変数名をセットすると、軸ラベルがage, weightと 表示されます。
# データのセット
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')
# グラフ表示
graph = preGraph("images/ch12_fig3.pdf")
r.plot(age, weight)
postGraph(graph)
Rのplot関数にタイトル等のオプションを指定する場合には、r関数を使って以下のように表示します。
# タイトル等グラフオプションを指定する場合には、r.plotではなくr関数を使用する
graph = preGraph("images/ch12_fig4.pdf")
# r('par(family="Japan1Ryumin")')
r("plot(age, weight, main='Test plot')")
postGraph(graph)
Rの得意とする集計処理を以下のようなデータファイルを読み込み、 クロス集計を作成してみましょう。
%%writefile data/ch12_data.txt
diabetes status
Type1 Poor
Type2 Improved
Type1 Excellent
Type1 Poor
Writing data/ch12_data.txt
データファイルの読み込みは、Rのread.table関数を使用します。 header=Tを付けて、1行目がヘッダであり、 1列目がdiabetes、2行目がstatusと名前付けします。
# データファイルの読み込み
d = r("d<-read.table('data/ch12_data.txt', header=T)") ; d
diabetes status 1 Type1 Poor 2 Type2 Improved 3 Type1 Excellent 4 Type1 Poor
クロス集計には、table関数を使用し、データセットdの変数、 diabetes, statusを指定します。
# クロス集計
r('table(d$diabetes, d$status)')
Excellent Improved Poor Type1 1 0 2 Type2 0 1 0