今、勉強中

chisq.test(d01, correct = FALSE)

x1<-rnorm(300) # ホワイトノイズ
x2<-rep(c(15,rep(0,49)),6) # スパイク異常値が6回生じる時系列で周期は50期
x3<-c(rep(0,150),rep(3,150)) # ステップ状の不連続な変化
x<-x1+x2+x3 # 全部合わせる
xt<-ts(as.numeric(x),frequency=50) # 必ずnumeric型で与えること！！！
# なお、frequency=50で50期ごとという周期の情報がインプットされる
plot(xt)

xt.stl<-stl(xt,s.window="periodic")
plot(xt.stl)

#http://www.cc.aoyama.ac.jp/~t41338/lecture/aoyama/stat2e/stat2e_top.html
#指導法データの csv ファイル（shidouhou.csv）　を以下のフォルダへダウンロード

setwd("d:/R/yasashi")
d01 <- read.csv("shidouhou.csv",head=T)
数学　<- d01$math
統計　<- d01$stat
table(数学,統計)
##結果
# 統計
# 数学 嫌い 好き
# 嫌い 10 4
# 好き 2 4

社会人のためのデータサイエンス演習を Rで

　まあ、なんと申しますかひと昔前だったら退職している年齢になってきてますので、お金の代わりに平日の時間にゆとりが生まれつつあります。
　ということで、お金のかからないアフターファイブなるものを追及しようとする今日この頃です。

　最近は、インターネットでいろいろな勉強ができるいわゆるＯＰＥＮな教育環境が整ってまいりました。私が最近受講を始めましたのは、gaccoというサイトの「ga063 社会人のためのデータサイエンス演習」です。

　この授業では、EXCELのデータ分析アドイン機能なんかをベースに課題を解いていくのですが、EXCELベースの分析作業というのは、ＧＵＩで試行錯誤していくうちにドロドロになってしまうどうも後味の悪い印象が強くて、Ｒの方がそういう意味ではスクリプトベースで再現性が保証されていて筋の良さを感じています。じゃあＲで分析してみたらどんなもんやいということをやってみたわけです。

　ちなみにRの操作能力は、業務等でバリバリ使っているわけでもありませんので、あくまで初学者プラスアルファ位という状況ですから、自分自身の備忘録という役割以上のものではありません。初学者が陥る前に進めないポイントについては、多少は役に立つことを書いているかもしれませんが、もしかすると、恐ろしく筋の悪い解決策を提示しているのかもしれません。
　非難・苦情・助言のたぐいは、コメントなり、メールなりで教えていただければ幸いです。

受講しないと（無料です）データとか取得できないし、課題の問題もわからないので、全然意味の分からない感じもいたしますけど、「わかるやつだけわかればいいという」心で書いています。



ｗｅｅｋ２分です。データは週ごとに別のフォルダに置いてあることとします。例えばweek2は、
D:\gacoo\data_sience\week2
というフォルダです。
------------------------------------------------------------------------------------------------------------

setwd("D:/gacco/data_sience/week2")
#数値の表示桁数を12ケタに


options(digits=12)

##　read.csv　で1行目が見出し情報の場合、head=T　とする。
## 数値が文字列で入ってくるので、文字列をファクター型にしないようにして
##　まず読み込む。
d02 <- read.csv ("dummydata_A.csv",head=T,stringsAsFactors=F)
#数値フィールドのカンマを取る
d02a <-sapply(d02[,c(10:21)],function(f) as.numeric(gsub(",","",f)))
#ファクター化したい項目だけ、変換する。
d02b <-sapply(d02[,c(2:8)],function(f) as.factor(f))
##一つのデータフレームに結合
d02c <- data.frame(d02[,1],d02b,d02[,9],d02a,stringsAsFactors=F)
#ヘダーの先頭が数値だったり、カッコが入っていると、きたない列名になるので
##　抜き出して編集してみる。
(d02c.name <- colnames(d02c))
#面倒なので、手で直したものを入れる
colnames(d02c) <- c("id","３大都市圏か否か","世帯人員","就業人員"
,"住居の所有関係","就業.非就業の別","年齢階級５歳階級","年齢階級６５歳未満か否か"
,"集計用乗率","年間収入","消費支出額_合計","食費","住居費","光熱水道費"
,"家具家事用品費","被服及び履物費","保健医療費","交通通信費","教育費","教養娯楽費"
,"その他の消費支出")

#d02c の構造
str(d02c)

##
##問題2-1　利用データdummydata_A.xlsx 消費支出のうち、食費の平均値、分散、標準偏差の
## 組み合わせとして正しいものを、次の (1)～(4)のうちから１つ選びなさい。
##
##　不偏分散というのは(n-1)で除（÷）したもの、標本分散は（n）で除したもの。
## この問題は前者で解答
##
mean(d02c$食費)
var(d02c$食費)
sd(d02c$食費)

##問題2-2 利用データ　ummydata_A.xlsx 階級間隔を10,000円として食費の度数分布表を
## 作成した時、20,001～30,000円以下の階級の度数として正しいものを選びなさい。

##　subsetでdata.frameから該当レコードを抜き出して、レコード数を表示
d02d <- subset(d02c, (　食費 >= 20001 & 食費 <= 30000 ))
nrow(d02d)

##問題2-3 階級間隔を10,000円とした食費のヒストグラムとして、最も近いものを選びなさい。

##　度数分布表をつくる
bunpu　<-table(cut(d02c$食費,
　　　　breaks = c(seq(0, to = 200000,by = 10000),max(d02c$食費))
　　　,dig.lab=10,include.lowest = T))
bunpu

#階級間隔を10,000円とした食費のヒストグラム
plot(bunpu)
##　線じゃ嫌なら
barplot(bunpu)
##ざっくりヒストグラムを描くなら
hist(d02c$食費,breaks=20)


#問題2-5 (3/3 点) 「その他の消費支出」を除く消費支出の中で最もばらつきの大きい項目を、
#　　次の(1)?(4)のうちから?つ選びなさい。（ばらつきの指標は標準偏差とする）
sd(d02c$食費)
sd(d02c$住居費)
sd(d02c$保健医療費)
sd(d02c$交通通信費)

R備忘録

怪しげなアマチュアスピリットにあふれるプログラマが、業務上Ｒに手をだして、まともに
使いこなせる訳が無いとは思うのだが、たぶんこの仕事EXCELではできなかったと思うよ。

まず、EXCEL-VBAでいうと、sub　と　functionで目的がちがうのだが、Ｒでは、
functionでぜんぶやるみたいだ。

↓あらら、"blockquote"タグを４つつなげたらこんな風になるんだ、、（追記）

MyFunc <- function(in_a,in_b) {
pdf(file=in_a,family="Japan1GothicBBB",width=8,height= 11.7)
plot(in_b,in_b)
dev.off()
# return (123)
}
MyFunc("d:/test.pdf",c(1,2,3,4,5,6,7))

切捨て：trunc()
剰余　：%%

Rでイールドカーブをグラフ描画してみる

Excelでは、グラフ印刷の場合、1行目、1列目が座標軸になる。
ビジネス系の場合、案外ｘ軸は非数値である場合がある。
そういう時のテクニック。
やや、冗長ですねえ、、もう少し、シンプルになると思います、、、
ポイントは、plotで「xaxt="n"」、 axis(1, at=<x軸の目盛り位置>, labels=＜その時のラベル＞, pos=0)

x <-read.csv("c:/test/test1.csv",header=T)
x.data<-data.frame(x[,2:14],row.names=x[,1])
x.data
lavel <-c(paste(c(1,3,6),"m",sep=""),paste(c(1:10),"y",sep=""))
plot(1:13,x.data[2,],lty=1,type="l",xaxt="n",xlab="イールド",ylab="利率（％）")
for(i in 1:13){
axis(1, at=i, labels=lavel[i], pos=0)
if (i>2){
lines(1:13,x.data[i,],lty=3,col=i)
}
}

> x.data
１ヶ月 ３ヶ月 ６ヶ月 １年 ２年 ３年
金利１ 0.100000000 0.30000000 0.50000000 0.70000000 0.90000000 1.10000000
金利２ 0.027877606 0.21092763 0.41189755 0.01992961 0.48301959 0.99195173
金利３ 0.072876887 0.04525560 0.35576517 0.08349966 0.40879922 0.95565290
金利４ 0.036667735 0.25436784 0.29419797 0.20294919 0.04106208 0.56739820
金利５ 0.005569511 0.25768794 0.36511799 0.35520448 0.44750703 0.61493886
金利６ 0.039964086 0.00021958 0.33824404 0.35769531 0.73470825 0.83909638
金利７ 0.085460670 0.04345718 0.08710587 0.62109500 0.05257187 0.29153159
金利８ 0.089940628 0.03554693 0.25882839 0.22154714 0.12280772 0.70943057
金利９ 0.072437949 0.24836105 0.28228370 0.46638298 0.68881066 0.08428973
金利１０ 0.009525149 0.01667500 0.07742818 0.30676718 0.37501316 0.83389939
４年 ５年 ６年 ７年 ８年 ９年
金利１ 1.3000000 1.5000000 1.7000000 1.90000000 2.1000000 2.30000000
金利２ 0.6803378 1.3745004 0.5135321 1.85668185 1.8031573 1.95940218
金利３ 0.2805778 0.6336938 0.1442276 0.34588735 0.5009916 0.46885884
金利４ 0.8763258 1.0030744 1.2785844 0.65323567 1.1720169 0.93684936
金利５ 0.2531363 0.9277034 0.1598404 0.16635548 2.0411710 0.02640765
金利６ 0.8336032 1.2527894 1.2140475 1.84418647 2.0969549 0.73039513
金利７ 0.9223572 0.5911564 0.9323400 0.40266169 0.2741270 0.11215907
金利８ 1.0514676 0.4819668 0.7789621 0.05189584 1.5780077 0.48349885
金利９ 0.4702907 0.1383999 0.6809487 0.64755723 0.1941831 0.43439839
金利１０ 1.1814694 0.7459063 0.2013630 1.74521164 1.3455539 2.22703265
１０年
金利１ 2.5000000
金利２ 1.1260667
金利３ 0.8738908
金利４ 0.4824346
金利５ 0.2104674
金利６ 0.7593693
金利７ 0.7010074
金利８ 0.8558189
金利９ 1.6905044
金利１０ 1.9910796

R言語の おやっ？

多分言語仕様が、小さなスクリプトを書きやすくしている。　C言語と同様に、中括弧をブロック構造の定義に
使っているが、改行を区切り文字として採用している分、慣れないと、おやって時がある。
セミコロンを使うこともできるが、1行に複数の命令を入れる時に使うぐらいかな。ある意味 BASICぽい感じもする。

##正しく判定される場合（１）
a <-c(FALSE)
if(a) {
print ("true")
print ("true")
} else {
print ("false")
}

##正しく判定される場合（２）
a <-c(FALSE)
if(a) { print ("true")
} else { print ("false")}

## 正しく判定される場合（３）
a<-c(FALSE)
if(a) { print ("true")} else { print ("false")}

##正しく判定される場合（４）
a <- TRUE
if(a) print ("true") else print ("false")

##elseがエラーになる場合（elseは、ブロック終わりと同じ行に無いとエラー）
a <-c(FALSE)
if(a) { print ("true")
}
else { print ("false")}

正規性の検定

x <- c(1:6)
shapiro.test(x)
wilcox.test(x)

x <-c(rnorm(1:100))
shapiro.test(x)
wilcox.test(x)

http://www.obihiro.ac.jp/~masuday/resources/stat/r_t-test02.html

Ｒグラフを並べて印刷する方法

#win.printの使い方
#最初にプリンタデバイスを開くので、印刷ダイアログがでてしまう。（連続処理には向かない）
#以下は、Ａ４サイズに印刷するサンプル
#
win.print(width=8, height = 11.7)
layout(rbind( c(1, 1, 2, 2),
c(3, 4, 5, 6),
c(7, 8, 9, 10),
c(11, 12, 13, 14)),heights=c(2,1,1,1))
layout.show(14)
dev.off()
#
# pdfへの印刷
#
pdf(file="c:/test.pdf",width =8 ,height =11.7)
layout(rbind( c(1, 1, 2, 2),
c(3, 4, 5, 6),
c(7, 8, 9, 10),
c(11, 12, 13, 14)),heights=c(2,1,1,1))
layout.show(14)
dev.off()
#　pdfのサンプルは
pdf(file="c:/test2.pdf",width =8 ,height =11.7)
layout(rbind( c(1, 1, 2, 2),
c(3, 4, 5, 6),
c(7, 8, 9, 10),
c(11, 12, 13, 14)),heights=c(2,1,1,1))
plot(sin)
plot(cos)
plot(tan)
plot(asin)
plot(acos)
plot(atan)
plot(sin)
plot(cos)
plot(tan)
plot(asin)
plot(acos)
plot(atan)
plot(acos)
plot(atan)
dev.off()

Ｒでレイアウト

layout(rbind(c(1, 1, 2, 2),
c(3, 4, 5,0),
c(6, 7, 8,0),
c(9, 10, 11,0)),heights=c(2,1,1,1), respect=TRUE)
layout.show(11)
layout(rbind(c(1, 1, 2, 2),
c(3, 4, 5,6),
c(7, 8, 9,10),
c(11, 12, 13,14)),heights=c(2,1,1,1))
layout.show(14)

Ｒによるグラフ作成の勉強

以下のＱＡを自分なりに回答してみる。近似線の多項式の当てはめ問題（模範解答はスプラインでしたけど）
任意のX軸と曲線との交点（Y値）の値の算出

金先生の04年8月～
間瀬先生の「R 基本統計関数マニュアル」p187～

par(mfrow=c(3, 2))
#質問、ブラックボックスのシステムで観測値を補完している。ｘ軸が12.5の時のy軸値はいくつなのか？
V <- c(84,80,76,68,59,50,44,38,35,33,32,34,37,41,49,63,73,81)
plot(V,type="b")
abline(v=12.5)

#模範解答　スプラインとして処理する（１）
x <- 1:18
V <- c(84,80,76,68,59,50,44,38,35,33,32,34,37,41,49,63,73,81)
plot(x, V,type="p")
ans <- spline(x, V, n=100)
lines(ans, col="red")
spline(x, V, n=100, xout=12.5)



#直交多項式を利用した例（２）
fm2=lm(V~poly(x,4))
plot(x,V)
lines(x,fitted(fm2),lty=2,col=4,lwd=2)
fm2

#lmを使った駄目な例（３）
x<-c(1:18)
V <- c(84,80,76,68,59,50,44,38,35,33,32,34,37,41,49,63,73,81)
fm3=lm(V~a+b*x+c*x^2+d*x^3+e*x^4)
#plot(x,V)
#lines(x,fitted(fm3),lty=2,col=3,lwd=2)
#fm3

#正しいlm（３）
fm3=lm(y~x+I(x^2)+I(x^3)+I(x^4))
plot(x,V)
lines(x,fitted(fm3),lty=2,col=1,lwd=2)
summary(fm3)
#nls関数（４）
x<-c(1:18)
V <- c(84,80,76,68,59,50,44,38,35,33,32,34,37,41,49,63,73,81)
fm4=nls(V~a+b*x+c*x^2+d*x^3+e*x^4,start=c(a=1,b=1,c=1,d=1,e=1),trace=T)
plot(x,V)
lines(x,fitted(fm4),lty=2,col=3,lwd=2)
fm4

Rで散布図（仮）

#csvファイルの読み込み、１列目をrow.namesに設定

d <- read.csv("f:/Program Files/R/R-2.9.1/data1.csv",header=T,row.names=1)
d
library(lattice)
xyplot(d[,1] ~ d[,3])
plot(d)

plot(d[,1]*100,d[,3]*100,type="n",xlab="住宅シェア（％）",ylab="総預金シェア（％）")
text(d[,1]*100,d[,3]*100,rownames(d))

plot(d[,3]*100,d[,1]*100,type="n",ylab="住宅シェア（％）",xlab="総預金シェア（％）")
text(d[,3]*100,d[,1]*100,rownames(d))

plot(d[,3]*100,d[,4],type="n",ylab="預金顧客数（人）",xlab="預金シェア（％）")
text(d[,3]*100,d[,4],rownames(d))

plot(d[,4],d[,2],type="n",xlab="預金顧客数（人）",ylab="住宅顧客数（人）")
text(d[,4],d[,2],rownames(d))

#df1　へ預金・貸金の契約者のみ移行する
df2 <- data.frame(cbind(d[,2],d[,4]))
#colnames(df2) <- c(colnames(d[2]),colnames(d[4]))
df2
cor(df2)
plot(df2)

df2.lm <-lm(X2~X1,data=df2)
summary(df2.lm)
plot(df2)
abline(df2.lm,lwd=2)
plot(df2.lm)


cor(d[,4],d[,2])
d.lm <-lm(d[,4]~d[,2],data=d)
summary(d.lm)
plot(d)

plot(d[,4],d[,2])
abline(df2.lm,lwd=2)


#df1　へ預金・貸金のシェアのみ移行する
df1 <- data.frame(cbind(rownames(d),d[,1],d[,3]),row.names=1)
colnames(df1) <- c(colnames(d[1]),colnames(d[3]))
df1

df1.d <-dist(df1)
round(df1.d)
df1.d
df1.hc<- hclust(df1.d)
summary(df1.hc)
plot(df1.hc)
df1.hc<- hclust(df1.d,method="ward")
summary(df1.hc)
plot(df1.hc)

#df1　へ預金契約者、預金・貸金のシェアのみ移行する
df1 <- data.frame(cbind(rownames(d),d[,1],d[,3]),row.names=1)
colnames(df1) <- c(colnames(d[1]),colnames(d[3]))
df1

df1.d <-dist(df1)
round(df1.d)
df1.d
df1.hc<- hclust(df1.d)
summary(df1.hc)
plot(df1.hc)
df1.hc<- hclust(df1.d,method="ward")
summary(df1.hc)
plot(df1.hc)

【解決】【統計環境Ｒ】散布図での漢字ラベル描画がうまくいかない

Windows版日本語グラフパッチ

Rをインストールしてあるディレクトリの
etcRdevgaを以下のように変更してください。
引用元：http://www.okada.jp.org/RWiki/

原因は不明ですけど、win2000の場合は、「Rdevgr」で指定するフォントを、
漢字フォントに変えないとダメみたいです。
引用元のファイルと中身を入れ替えればｏｋでした。

（メモ）
・Ｒの2.51版で比較したところ、設定ファイル（Rdevgr）の項目は同じでした。
・フォントは「MS Mincho」から「MS Gothic」に置き換えしました。
　（これは趣味の問題でしょう）

【統計環境Ｒ】散布図での漢字ラベル描画

Ｒという統計環境をいま勉強してます。会社で動作が不安定だった
散布図の漢字ラベルですけど、自宅だと普通に上手く表示できる。
（会社では漢字が白い四角）
違いは、ＯＳの違い位と思うんだが、WIN2000（会社）でダメ
でＸＰ（自宅）でＯＫというだけなのかな？
もう少し調べてます。

#Book1.txt　の読み込み
data1<-read.table("Book1.txt",header=TRUE,row.names=1)
#表示
data1
##表示結果
# 数値１ 数値２
#東京 1 1.5
#大阪 2 2.5
#名古屋 3 3.5
#福岡 4 4.5
#札幌 5 5.5

#デフォルトの散布図
plot(data1[,1],data1[,2])
#漢字ラベルでの散布図描画
data1.label=c("東京","大阪","名古屋","福岡","札幌"）
plot(data1[,1],data1[,2],type="n")
text(data1[,1],data1[,2],data1.label)

#もし、漢字名でなくシーケンスを表示させるなら、
plot(data1[,1],data1[,2],type="n")
text(data1[,1],data1[,2])