統計学2020 Lecture 10

Rを用いた演習 (1) Rのインストールと基本文法

• 今回からRによる演習を行います．
  
• フォルダーをダウンロードしzipファイルを解凍します（zipフォルダーの中で作業をしないでください．
  
• Rを利用する際には，Lecture10_ToStudentフォルダーの「Rdata」というファイルをダブルクリックして開始してください．拡張子を見えなくしている場合には「R」と書かれたファイルがそれにあたります．
  
• 7月15日もリモートで質問を受け付けます (13:00-14:00とします）.

1 R事始め

1.1 Rのインストール

Windowsの場合，Rの下記サイトで“Download R for Windows”を選択し，“Base”の“Download R 4.x.x for Windows”を選択してインストーラーファイルをダウンロードすれば，自宅などのPCに簡単にインストールできます． (Macも同様です)

https://cran.r-project.org/

図1(a) Rのダウンロードサイト

図1(b) Baseを選択

図1(c) インストーラーをダウンロードし実行

1.2 Rstudioのインストール

Rを使うための大変便利なRstudioというソフトがあります．ごく稀にインストールがうまくできないことがありますので，今回はRをそのまま使うことを前提で説明を進めます．ただ興味のある方は，試しに使っていただければと思います．

https://rstudio.com/products/rstudio/download/#download

図2. Rstudioのサイト

2 Rを用いた簡単な計算

• デスクトップのRのアイコンをダブルクイックしてRを立ち上げる
  
• ためしに以下の計算をしてみる

入力の仕方は2通り：

• 直接Rのコンソール画面に入力（コードを打ち間違えたときに，いちからやり直しのが面倒．打ち間違えて入力し直すときにはEscキーを押していったん無効にしてから再度入力してください）

図3. 直接入力する場合

• R付属のエディターにコードを入力し，実行したいコードを選択して，「コントロール+R」で実行．

図4(a). Rのエディターを利用する場合

図4(b). 実行したいコードを選択して，「コントロール+R」で実行

2.1 基本演算

1+3*2-4/2   #四則演算

[1] 5

11^2        #べき乗

[1] 121

1/sqrt(4)   #平方根(square root)の逆数

[1] 0.5

4^{-1/2}    #上記と同じ計算

[1] 0.5

log(9)      #自然対数

[1] 2.197225

exp(6)      #指数

[1] 403.4288

sin(pi/4)   #三角関数

[1] 0.7071068

sqrt(2)/2   #上記の確認

[1] 0.7071068

a <- 3      #新しいオブジェクト（定数）の定義
a           #定義した定数の確認

[1] 3

a^2 + log(a) #それを用いた計算例

[1] 10.09861

2.2 ベクトル（配列）の扱い

a <- c(1,3,5,7,9)   #新たにベクトルをオブジェクトとして定義
a                   #定義したベクトルの確認

[1] 1 3 5 7 9

sum(a)              #合計

[1] 25

mean(a)             #平均

[1] 5

a+10                #ベクトルのすべての要素に10を加える

[1] 11 13 15 17 19

5*a                 #ベクトルのすべての要素を10倍する

[1]  5 15 25 35 45

a^2                 #ベクトルのすべての要素を2乗する

[1]  1  9 25 49 81

log(a)              #ベクトルのすべての要素の対数を取る

[1] 0.000000 1.098612 1.609438 1.945910 2.197225

a[3]          #ベクトルaの3番目の要素を取り出す

[1] 5

a[c(1,2,3)]   #ベクトルaの1,3,5番目の要素を取り出す

[1] 1 3 5

a[-c(2,4)]    #ベクトルaの2,4番目以外の要素を取り出す（上記と同じ）

[1] 1 5 9

b <- seq(10,30,5) #同じ長さのベクトルbを新たに定義("seq(x,y,z)"はxからyまでz刻みで)
a+b  #要素どおしの足し算

[1] 11 18 25 32 39

a*b  #要素どおしの掛け算

[1]  10  45 100 175 270

b/a  #要素どおしの割り算

[1] 10.000000  5.000000  4.000000  3.571429  3.333333

2.3 繰り返し計算

1から10までの足し算を行う計算例．

ss <- 0
 for(i in 1:10){
  ss <- ss + i  
}
ss

[1] 55

a<- seq(1,10,1)
sum(a)

[1] 55

3 外部データの読み込みと可視化

気温の年変化を図で表してみよう．東京の気温と降水量のデータを“Data_Tokyo.csv”というファイルに保存してあるとします．

3.1 データの読み込み

Data_Tokyo <- read.csv("Data_Tokyo.csv", header=T) #CSVファイルの読み込み
head(Data_Tokyo, 10) #読み込みの確認（最初の10行）

   Year Month Temperature Rainfall
1  1900     1         1.5     69.2
2  1900     2         2.6     31.3
3  1900     3         5.2     69.3
4  1900     4        12.3    168.6
5  1900     5        17.5    127.3
6  1900     6        20.4     89.4
7  1900     7        23.8    122.1
8  1900     8        26.6     66.3
9  1900     9        22.6    173.9
10 1900    10        15.4    113.6

tail(Data_Tokyo, 10) #読み込みの確認（最後の10行）

     Year Month Temperature Rainfall
1431 2019     3          NA       NA
1432 2019     4          NA       NA
1433 2019     5          NA       NA
1434 2019     6          NA       NA
1435 2019     7          NA       NA
1436 2019     8          NA       NA
1437 2019     9          NA       NA
1438 2019    10          NA       NA
1439 2019    11          NA       NA
1440 2019    12          NA       NA

names(Data_Tokyo)    #Data_Tokyoのオブジェクトの構成要素を確認

[1] "Year"        "Month"       "Temperature" "Rainfall"   

dim(Data_Tokyo)      #Data_Tokyoの配列次元の確認1392×4の行列であることが分かる

[1] 1440    4

データはData_Tokyoという名前のオブジェクトに保存されている．ここから気温や降水量を取り出したい時には，Data_Tokyo$Temperature'やData_Tokyo$Rainfall’として取り出せばよい．あるいはData_Tokyo[,3]とすれば気温を取り出せる．また1月の気温だけを取り出す場合には，下記のようにすればOK.

Data_Tokyo$Temperature[Data_Tokyo$Month==1]

 [1] 1.5 4.2 1.7 3.8 1.5 4.3 1.6 3.8 3.2 2.5 4.4 2.8 2.3 1.8 3.5 3.2 4.8 1.4 1.0
[20] 2.9 3.3 3.6 1.3 1.5 2.7 2.4 2.2 2.7 3.6 1.7 2.7 3.4 4.6 2.7 1.1 3.2 0.9 4.0
[39] 2.2 2.0 2.3 4.6 3.0 1.8 2.9 0.9 3.4 3.8 3.5 4.6 5.2 3.5 4.1 2.9 5.1 3.2 3.6
[58] 4.3 4.0 3.4 3.9 2.6 3.2 1.1 5.5 3.3 3.4 3.2 3.4 4.8 3.6 4.3 7.1 6.2 3.2 4.0
[77] 3.7 2.5 4.7 5.7 4.2 2.4 3.8 4.8 2.2 3.1 2.8 4.8 5.9 6.8
 [ reached getOption("max.print") -- omitted 30 entries ]

毎回ドルマークを付けて取り出すのが面倒な場合は，以下のように最初にattachという関数を利用してData_Tokyoの名前を覚えさせてから実行してもよい（あまり乱用しない方がよいやり方ではあるが）．

attach(Data_Tokyo)
Temperature[Month==1]

 [1] 1.5 4.2 1.7 3.8 1.5 4.3 1.6 3.8 3.2 2.5 4.4 2.8 2.3 1.8 3.5 3.2 4.8 1.4 1.0
[20] 2.9 3.3 3.6 1.3 1.5 2.7 2.4 2.2 2.7 3.6 1.7 2.7 3.4 4.6 2.7 1.1 3.2 0.9 4.0
[39] 2.2 2.0 2.3 4.6 3.0 1.8 2.9 0.9 3.4 3.8 3.5 4.6 5.2 3.5 4.1 2.9 5.1 3.2 3.6
[58] 4.3 4.0 3.4 3.9 2.6 3.2 1.1 5.5 3.3 3.4 3.2 3.4 4.8 3.6 4.3 7.1 6.2 3.2 4.0
[77] 3.7 2.5 4.7 5.7 4.2 2.4 3.8 4.8 2.2 3.1 2.8 4.8 5.9 6.8
 [ reached getOption("max.print") -- omitted 30 entries ]

3.2 気温データの図示（標準的な方法）

#横軸に年，縦軸に気温，それを"line"で結ぶ
plot(Year, Temperature, type="l", main="東京の気温変遷", xlab="年", ylab="気温(°)")  

月毎にプロットすると変化の様子が分かりやすい．

par(mfrow=c(3,4), mar=c(2,2,2,2)) #図を3×4で表示

MEAN <- numeric(12)  #月別の平均気温を計算して保存するためのオブジェクト
for(i in 1:12){
 MEAN[i] <- mean(Temperature[Month==i], na.rm=T) #月別の平均気温の計算
 plot(Year[Month==i], Temperature[Month==i], type="l") #月別の気温変化プロット
 abline(MEAN[i], 0, col="red", lty=2)  #切片MEAN[i],傾き0の赤い直線(鎖線)
 mtext(paste("Month", i), side=3, cex=0.7) #月を記入
}

3.3 気温データの図示（自動的に綺麗な図を書く方法）

ggplot2というパッケージを使うと図が大変きれいに書けます．ちょっと文法を覚えるのが面倒ですが．

始めて使う場合にはまずパッケージをインストールします．

install.packages("ggplot2")   #インストール（1回だけ実行すればOK）
library(ggplot2)              #パッケージの読み込み（Rを立ち上げるごとに毎回実行）

先ほどの月別気温の変化をカラフルに書いてみます．

tmp <- ggplot(Data_Tokyo, aes(Year, Temperature)) + geom_line()
tmp 

df <- data.frame(Month,MEAN)
tmp <- ggplot(Data_Tokyo, aes(Year, Temperature, col=factor(Month))) + 
   geom_line() + facet_wrap(.~Month, ncol=4) + 
   geom_hline(data=df, aes(yintercept=MEAN), linetype="dashed", color="red")
tmp

縦軸のスケールを変えが方が年変化の様子が分かりやすく，また凡例もこの場合は不要なので削除．

tmp <- ggplot(Data_Tokyo, aes(Year, Temperature,col=factor(Month))) + 
   geom_line(show.legend=FALSE) + facet_wrap(.~Month, ncol=4, scales="free") +
   geom_hline(data=df, aes(yintercept=MEAN), linetype="dashed", color="red")
tmp

3.4 降水量データの図示（標準的な方法）

plot(Year, Rainfall, type="l", main="東京の降水量変遷", xlab="年", ylab="気温(°)")  

par(mfrow=c(3,4), mar=c(2,2,2,2)) 

MEAN <- numeric(12)  
for(i in 1:12){
 MEAN[i] <- mean(Rainfall[Month==i], na.rm=T) 
 plot(Year[Month==i], Rainfall[Month==i], type="l") 
 abline(MEAN[i], 0, col="red", lty=2) 
 mtext(paste("Month", i), side=3, cex=0.7) 
}

df <- data.frame(Month,MEAN)

tmp <- ggplot(Data_Tokyo, aes(Year, Rainfall, col=factor(Month))) + 
   geom_line(show.legend=FALSE) + facet_wrap(.~Month, ncol=4, scales="free") +
   geom_hline(data=df, aes(yintercept=MEAN), linetype="dashed", color="red")
tmp

4 提出課題（提出日は後日お伝えします）

気象庁のHPから自分の出身地のデータ（月別気温と降水量）をダウンロードし，上記のように気温と降水量の時系列プロットを作成する事．（標準的な図あるいはggplotのどちらでもよい）

4.1 気象庁HPを訪問

アドレス http://www.data.jma.go.jp/gmd/risk/obsdl/index.php

4.2 都道府県を選択

4.3 さらに場所を選択

4.4 月別の気温を指示

4.5 月別の降水量を指示

4.6 期間を設定(1900年1月から2019年12月)し，CSVファイルをダウンロード

4.7 データファイルの作成

Data_Tokyo.csvをコピーし，名前を付け替える．Data_Kyoto.csvなど．そのファイルをエクセルで開けて，ダウンロードした気温と降水量データをコピペする．あとは，上記と同じように作業をすれば図が作成できます．トライしてみて下さい．提出物は計算に利用したコードと図です．