いろいろと名言があった。
人間はだれでも、人間としての存在の完全なかたちを備えている。
世界でいちばん高い玉座の上にあがったとしても、われわれはやはり、自分のお尻の上に座るしかない。
など。
日経平均とドル円と売買代金の分析２- R言語で平均値の標準誤差(Standard Error)と信頼区間(Conficence Interval)を求める。
の続きです。
今回は、Michael J. CrawleyのStatistics: An Introduction Using Rを参考にして、日経平均、ドル円、売買代金の月次変化率の平均値の標準誤差と信頼区間を求めてみようと思います。
Statistics: An Introduction Using R
- 作者:Michael J. Crawley
- 出版社/メーカー: Wiley
- 発売日: 2014/11/24
- メディア: ペーパーバック
上の画像のような本です。データのCSVファイルをウェブからダウンロードできるので実際にR言語を使って自習できます。
variance(分散)はvar関数で取得できます。
日経平均の前月比のvarianceは0.001445938、ドル円の前月比のvarianceは0.0003694305、売買代金の前月比のvarianceは0.01295124です。売買代金のvarianceが一番大きいですね。
これらのvarianceが統計的に有意に違うのかどうか、var.test関数で調べてみましょう。
日経平均とドル円のvarianceを比較します。
p-value = 5.962e-07と0.05よりも小さい値なので、日経平均の前月比のvarianceとドル円の前月比のvarianceは有意に違いがあることがわかります。ということは、売買代金の前月比のvarianceは他の二つとはもっと離れていますから、それぞれのvarianceは有意に違うのですね。
それでは標準誤差(standard error)を計算しましょう。
標準誤差 = SQRT(variance / n)です。nはデータの数です。
length(df$ChgNikkei) - 1 と1を引いているには、NAが一つあるからですね。
それぞれの平均値も計算しましょう。mean関数です。
Chgのデータ数(NAを除いた数) を数えましょう。length関数とis.na関数とsum関数を使うと計算できます。
データの数は59個ですね。
標準誤差と平均値とデータの個数を使ってこのように記述することができます。
日経平均の月次変化率の平均値は、1.005526 ± 0.004950497(1 s.e., n = 59)
ドル円の月次変化率の平均値は、0.9985149 ± 0.002502306(1 s.e., n = 59)
売買代金の月次変化率の平均値は、1.006599 ± 0.01481595(1 s.e., n = 59)
なんか学術論文っぽいですね。笑
プラスマイナスだと、いまいちイメージがつかめないので、実際に範囲を計算してみます。
まず、マトリックスを作成しました。
これで三つを一度に計算できます。
日経平均は平均値 - 1標準偏差の値が1.0005754と1よりも大きいですので、分析期間中は上昇傾向だったのかもしれないですね。ドル円と売買代金は 1をまたいでいますから分析期間では上昇傾向でも下降傾向でもなかったということですね。
信頼区間も計算します。
n = 59ですから自由度は59 -1 = 58 です。95%の信頼区間は、qt関数を使って、qt(0.025, 58)、qt(0.975, 58)を標準誤差にかけた値になります。
95%の信頼区間(n = 59)では、日経平均の月次変化率の平均値は、0.9956164から1.015435の範囲になり、1.0をまたぎますので、上昇傾向、下降傾向、どちらともいえないですね。
マトリックスでは、mtx[ , "AVG"]のように列名で指定することもできます。
今回は以上です。
日経平均とドル円と売買代金の分析１- 売買代金が一番変動が大きく、ドル円が一番変動が小さい。
今回は、日経平均とドル円と売買代金のデータを分析してみようと思います。
日経新聞電子版からデータをCSVファイルで取得しました。
こういうデータです。Nikkeiが日経平均、Yenがドル円、Daikinが一日当りの売買代金平均(10億円)です。2014年12月から2019年11月までのデータがあります。
これをR言語で分析したいと思います。
read.csv関数でCSVファイルを読み込み、head関数ではじめの6行を表示したのが上の画像です。
今回はこのままのデータではなくて、前月比を分析したいと思います。
前月比を作らないといけないので、このデータフレームにそれぞれの前月の値を追加します。
一番初めの2014年12月は当然、前月のデータが無いのでNAになります。ので、c(NA ~~~)とします。そして、最後のデータだけ削除したします。しdf$Nikkei[-length(df$Nikkei)]をくわえると、上の画像のようになります。
そして前月比を計算しましょう。
summary関数で各変化率のサマリを出力してみましょう。
売買代金が一番変動が大きいようです。それぞれの変動係数を計算してみます。
sd関数で標準偏差を計算して、mean関数で平均値を計算し、標準偏差 / 平均値で変動係数を算出しています。やはり、売買代金の変動係数が一番大きいですね。ドル円が一番変動が小さいです。
今回は以上です。
国立科学博物館の標本データベースの日本の鳥類の標本の分析５ - ウグイスとスズメを比べると、スズメは北海道が多い。１～３月はスズメが多く、１０～１２月はウグイスが多い。
の続きです。
今回はウグイスとスズメについて調べてみます。
まずは、ウグイス、スズメだけの作業用のデータフレームを作成します。
ウグイスは124、スズメは121の標本があります。
それでは、都道府県別の標本数を見てみます。table関数を使います。
う～ん、どうなんでしょうか。。正直この表を眺めてもよくわからないですね。
barplot関数を使って棒グラフを描きます。
どうでしょうか？
都道府県の数が多いのでもう少し集約しましょう。
これでもう一度、都道府県別に集約してみます。
これだと、ウグイスは鹿児島が多く、スズメは北海道が多いとわかりますね。
カイ二乗検定をしてみます。chisq.test関数です。
p値が0.005635ですから、ウグイスとスズメでは違いがあるということですね。
上の表のどの箇所が有意に違うのか見てみます。
調整済み残差を見てみます。kaiというカイ二乗検定の結果の中に、stdresという名前で格納されていますので、kai$stdresで呼び出せます。
調整済み残差で見ると、絶対値で2以上のところが統計的に有意な違いがあるところです。その他と北海道ですね。スズメはウグイスと比較すると、北海道で採集された標本が有意に多いということですね。
年別ではどうでしょうか？taable関数です。
この表を見ただけではよくわからないですね。barplot関数を使って棒グラフにしてみます。
棒グラフにしてみました。スズメは1930年代、ウグイスは1970年か80年代に多いような気がします。年代で区切ってみましょう。
round関数で引数を-1にすると一桁上の数で丸めてくれます。-4.9を引いてから一桁上で丸めれば10年ごとの年代別になります。すこしわかりやすくなりました。ウグイスはどの年代にも標本がありますが、スズメは1890年代、1940年代、1980年代は標本がありません。
月別にも集計してみます。
これも一見しただけではよくわからないですね。棒グラフを描きます。
１、２、３月はスズメが多くて、１０、１１、１２月はウグイスが多いようですね。
12か月を１～３、４～９、１０～１２の３つに分けて集計してみましょう。
こうして３つのカテゴリーにしました。
棒グラフにしてみます。
4月から9月は同じくらいの差で、１、２、３月はスズメ、１０、１１、１２月はウグイスが多いとはっきりわかりますね。この差は有意な差なのか、カイ二乗検定をしてます。
p値が0.0006716です。0.05よりも有意な差があるということですね。
今回は以上です。
国立科学博物館の標本データベースの日本の鳥類の標本の分析４ - メジロは東京、ヒヨドリは東北の標本が多く、2,3,4月はヒヨドリが多く、1973年はヒヨドリがとっても多かった。
の続きです。
今回はメジロとヒヨドリでは、標本の年、月、都道府県に差があるか調べてみます。
まずは、メジロとヒヨドリだけの作業用のデータフレームを作ります。
Nameには、アオアシシギなど不要のファクタ水準があるので、削除します。
as.character関数で文字列にして、as.factorでファクタに戻しています。
では、まずは、都道府県を見てみましょう。
愛知県とか愛媛県とかメジロもヒヨドリも0の都道府県が何個かありますね。整理しましょう。
すこしすっきりしました。眺めると、ヒヨドリは岩手、宮城、秋田、青森など東北が多くて、メジロは東京が多いですね。
それぞれの東北で採取された標本数を調べます。
まずは、東北と東北以外というファクタを作成します。
岩手県、宮城県、山形県、秋田県、青森県を東北にして、それ以外を東北以外にしました。これでクロス表を作ります。
クロス表を作成したら、chisq.test関数でカイ二乗検定をします。
p値が9.184e-05ですから東北か東北でないかはヒヨドリ、メジロの標本数に関係があります。ヒヨドリは東北で採取された標本が多いということですね。
同じ要領で東京ファクターを作成してみましょう。
このように、東京以外と東京都にわけました。
それではやってみます。
p値は0.01344で0.05よりも小さいので、メジロとヒヨドリでは東京か東京以外かという標本数に有意な違いがある、ということですね。
今度は年に違いがあるかどうかを調べます。
どうでしょうか。。。これは年の数が多すぎてよくわからないですね。。
barplot関数でグラフにしてみます。
1973年、ヒヨドリは36羽も標本になっています。突出して多いですね。何かあったんでしょうか？
月別はどうでしょうか？
なんとなくヒヨドリとメジロでは採集された月に違いがある感じがします。
barplot関数で視覚化します。
２、３、４月はヒヨドリが多いですが、それ以外の月はメジロが多いですね。
２、３、４月とそれ以外の月で二分してクロス表にしてみましょう。
これでクロス表を作成してカイ二乗検定をします。
p値が6.822e-06と0.05より小さいのこのクロス表には有意な偏りがあります。ヒヨドリはメジロと比べると2,3,4月に採集された数が有意に多いことがわかります。
今回は以上です。
国立科学博物館の標本データベースの日本の鳥類の標本の分析３ - Mokuのデータが怪しいことに気が付いた。
の続きです。
今回は、Moku(目), Zoku(属), Shu(種)について調べます。
まず、Mokuは何種類あるか、標本数の多いMokuは何か調べてみましょう。
levels関数とlevels関数でMokuの種類数を調べます。
Mokuは214種類ですね。どのMokuが一番多いか、summary関数で見てみます。
sort関数で小さい順に並び替え、rev関数で大きい順にして、head関数で上位6を表示しています。PasseriformesというMokuが一番多いです。
Zokuも同じように調べます。まずは、何種類のZokuがあるでしょうか？
232種類です。
どの属が一番多いでしょうか？
EmberizaというZokuが一番多くて、557個ありました。
同じように、Shu(種)を調べましょう。
Shuの数は382種類、japonicusというのが212で一番多いです。おそらくこれがメジロなのでしょう。確認します。
idx <- df$Shu == "japonics"という行で論理ベクトルをつくり、これでdfの中からShuがjapinicusのだけを検索しています。head関数で始めの6行だけ表示しています。やっぱりメジロですね。メジロのMoku,ZokuはZosteroposというので、Mokuの一番多いPasseriformesでは無いのですね。
二番目に多いShuのamaurotisは何という鳥でしょうか？
amaurotisはヒヨドリのことでしたね。ヒヨドリもMokuは一番多いPasseriformesではないですね。
それでは、Passeriformesにはどんな鳥がいるか見てみましょう。
あれ？メジロ、ヒヨドリが入っています。おかしいですね。。
もとのCSVファイルをもう一度みてみましょう。
こんなふうにメジロでもZosteriposになっているのもあれば、Passerifermosになっているのもありました。なんでですかね。。。Mokuはちょっと信用できないってことですね。Zokuはどうでしょうか？一番数の多かったEmberizaを見てみます。
Emberizaに属するのはアオジ、ホオジロ、カシラダカなどとわかりました。
今回は以上です。
国立科学博物館の標本データベースの日本の鳥類の標本の分析２ - 月と日付のデータを整える。８月に採取された標本が一番少ない。
の続きです。今日は月と日付のデータを整えます。
Monthを確認します。
空白、**, 20は明らかにおかしいですね。NAにします。空白は1番目、**は2番目、20は8番目です。
こうして１から１２までにしました。これを文字列型にしてから数値型にします。
summaryの結果が度数表ではなく、平均値などになっています。数値型になりました。
Dayも同じようにして整えます。
Dayにも空白と**がありますね。これをNAにします。
はい、空白と**がなくなりました。これを同じように文字列に変換して数値型に変換します。
はい。できました。
それではMonthをヒストグラムにしてどの月に採集された標本が多いか調べます。
自作したdosuu関数を使います。
10月、11月、7月が多いですね。8月が一番少ないです。8月は夏休みだからみんな採集に出かけるかと思っていましたがそうじゃないんでしょうか？
Dayもヒストグラムで見てみましょう。自作したdosuu関数を使います。
バラバラって感じですね。
今回は以上です。
