Hatena::ブログ(Diary)

goinger的日記

2010年09月07日

[R] 多変量時系列解析メモ

わりかしRとかで一変量(というべきか..), ARIMA, GARCHあたりの実行手法を解説してる本はあるんですけど、多変量時系列(VAR, VARIMA)あたり解説してる本は無いのでメモ。

library(vars)

US.var <- VAR(cbind(GNP, M1), p = 3, type = "trend")
coef(US.var)

acf(resid(US.var)[, 1])
acf(resid(US.var)[, 2])

US.pred <- predict(US.var, n.ahead = 4)
US.pred

plot(US.pred)

GNP.pred <- ts(US.pred$fcst$GNP[,1], st = 1988, fr = 4)

M1.pred <- ts(US.pred$fcst$M1[,1], st = 1988, fr = 4)

ts.plot(cbind(window(GNP, start = 1981), GNP.pred), lty = 1:2)
ts.plot(cbind(window(M1, start = 1981), M1.pred), lty = 1:2)

まあ必要な人はあんまいないのではないかと

2010年08月01日

[]多変量解析手法の簡易メモなど

Rによるデータサイエンスって本のデータ解析メモですね.

この本はわりかし一般的な手法の網羅的な解説 + 参考文献豊富で良い感じ.


主成分分析(Principal Component Analysis)

  • 目的
    • 多変量データを少ない変数で表現できるようにする
    • 通常は2~3変数に縮約する場合が多い。
    • 分散を最大化する手法
  • 多変量解析としては最も有名な手法の一つ
  • 分散共分散行列の固有値問題とみなす
    • 最も固有値が大きいのが第一主成分、次に固有値が大きいのが第ニ主成分... と求める
    • 主成分の解析には寄与率、累積寄与率に関する情報を利用.
  • 主成分得点、予測、biplot...
  • R
    • princompやprcompなどが利用可能
    • princompを普通利用する
    • prcompはデータのスケールなど利用可能

サンプル

require(graphics)

## The variances of the variables in the
## USArrests data vary by orders of magnitude, so scaling is appropriate
(pc.cr <- princomp(USArrests))  # inappropriate
princomp(USArrests, cor = TRUE) # =^= prcomp(USArrests, scale=TRUE)
## Similar, but different:
## The standard deviations differ by a factor of sqrt(49/50)

summary(pc.cr <- princomp(USArrests, cor = TRUE))
loadings(pc.cr)  ## note that blank entries are small but not zero
plot(pc.cr) # shows a screeplot.
biplot(pc.cr)

## Formula interface
princomp(~ ., data = USArrests, cor = TRUE)
# NA-handling
USArrests[1, 2] <- NA
pc.cr <- princomp(~ Murder + Assault + UrbanPop,
                  data = USArrests, na.action=na.exclude, cor = TRUE)
pc.cr$scores

因子分析

  • 目的
    • 主に心理学社会学などで外的基準がない量的データから共通因子を見つけ出す探索的データ解析時に利用
    • 変数間の相関関係から共通因子を求める
  • 歴史
    • 1904 Spearmanによって提唱された
  • 利用の際の注意点
    • 多義的解釈が可能なので、自分に都合の良い解釈が可能なので注意が必要
    • 客観的意味付けをできるように使う事が重要
    • 多義性が少ない主成分分析、対応分析などを兼用すると吉
  • アルゴリズム
    • 主因子法、最尤法などなど
      • 主因子法: 安定した結果が得られる
      • 最尤法 : データが正規分布に従うときに利用すると吉。但し初心者向きではない
  • 因子の回転
    • 解釈の便宜のため,高い相関をもつ項目を共通因子として空間上の軸を決める操作を行う。これが因子の回転である。
    • 種類
      • 直行回転、斜交回転
      • 直行回転: バリマックス(よく使われる)、バイコーティマックス、コーティマックス、エクィマックス
      • 斜交回転: プロマックス(よく使われる)、コバリミン、バイコーティミン、コーティミン

  • R
    • パッケージ: stats
    • 関数: factanal

サンプル

# A little demonstration, v2 is just v1 with noise,
# and same for v4 vs. v3 and v6 vs. v5
# Last four cases are there to add noise
# and introduce a positive manifold (g factor)
v1 <- c(1,1,1,1,1,1,1,1,1,1,3,3,3,3,3,4,5,6)
v2 <- c(1,2,1,1,1,1,2,1,2,1,3,4,3,3,3,4,6,5)
v3 <- c(3,3,3,3,3,1,1,1,1,1,1,1,1,1,1,5,4,6)
v4 <- c(3,3,4,3,3,1,1,2,1,1,1,1,2,1,1,5,6,4)
v5 <- c(1,1,1,1,1,3,3,3,3,3,1,1,1,1,1,6,4,5)
v6 <- c(1,1,1,2,1,3,3,3,4,3,1,1,1,2,1,6,5,4)
m1 <- cbind(v1,v2,v3,v4,v5,v6)
cor(m1)
factanal(m1, factors=3) # varimax is the default
factanal(m1, factors=3, rotation="promax")
# The following shows the g factor as PC1
prcomp(m1)

## formula interface
factanal(~v1+v2+v3+v4+v5+v6, factors = 3,
         scores = "Bartlett")$scores

## a realistic example from Bartholomew (1987, pp. 61-65)
utils::example(ability.cov)

対応分析

  • 対応分析とは
    • 頻度データ、質的データの個体と変数との関連性、パターン分析を行う手法
    • コレスポンデンス分析とも呼ばれる
    • 数量化III類と似ている

サンプル

farms.mca <- mca(farms, abbrev=TRUE)
farms.mca
plot(farms.mca)

多次元尺度法(MDS: Multi-Dimensional Scaling)

  • 目標
    • データの個体間の類似度や距離を求めてそれを2~3次元にプロットしてデータの構造やパターン形成などを把握する手法
  • 分類
    • 計量、非計量の二種類
  • 計量的MDS
    • 解析の流れ
      • 距離を求める
      • 座標値を求める
      • 2~3次元上で個体を配置する(散布図作成)
      • 信頼性などの考察
  • Rでの求め方
    • cmdscaleを利用

サンプル(ヨーロッパの都市の距離データ)

require(graphics)

loc <- cmdscale(eurodist)
x <- loc[,1]
y <- -loc[,2]
plot(x, y, type="n", xlab="", ylab="", main="cmdscale(eurodist)")
text(x, y, rownames(loc), cex=0.8)

cmdsE <- cmdscale(eurodist, k=20, add = TRUE, eig = TRUE, x.ret = TRUE)
utils::str(cmdsE)

  • 非計量的MDS
    • 計量的MDSは直接距離などを求められる事を前提てしているが、心理データなどの親近性データは距離の性質を満たさない
    • 距離の性質を満たさない類似性データも利用可能にしたのが非計量的MDS
  • R
    • パッケージ: MASS, mlbench, e1071, vegenなど
    • 関数: isoMDS, sammon, metaMDS

isoMDS

swiss.x <- as.matrix(swiss[, -1])
swiss.dist <- dist(swiss.x)
swiss.mds <- isoMDS(swiss.dist)
plot(swiss.mds$points, type = "n")
text(swiss.mds$points, labels = as.character(1:nrow(swiss.x)))
swiss.sh <- Shepard(swiss.dist, swiss.mds$points)
plot(swiss.sh, pch = ".")
lines(swiss.sh$x, swiss.sh$yf, type = "S")

クラスタ分析

  • 大きく分けて大きく三種類ある
    • 階層、非階層(グループ数指定)、モデルに基づく手法

階層的クラスタリング

  • 注意事項
    • 個体数が大きいと計算量が膨大になるので大規模データには不向き
    • 小規模データでも適用は計画的に
  • 解析の流れ
    • データから距離(or 類似度)を求める
    • クラスタ分析手法の適用
    • コーフェン行列を求める
    • コーフェン行列から樹形図を作る
    • 結果についての検討を行う
  • 手法
    • 最近隣法, 最遠隣法, 群平均法, メディアン法, 重心法, ウォード法

サンプル

require(graphics)

hc <- hclust(dist(USArrests), "ave")
plot(hc)
plot(hc, hang = -1)

## Do the same with centroid clustering and squared Euclidean distance,
## cut the tree into ten clusters and reconstruct the upper part of the
## tree from the cluster centers.
hc <- hclust(dist(USArrests)^2, "cen")
memb <- cutree(hc, k = 10)
cent <- NULL
for(k in 1:10){
  cent <- rbind(cent, colMeans(USArrests[memb == k, , drop = FALSE]))
}
hc1 <- hclust(dist(cent)^2, method = "cen", members = table(memb))
opar <- par(mfrow = c(1, 2))
plot(hc,  labels = FALSE, hang = -1, main = "Original Tree")
plot(hc1, labels = FALSE, hang = -1, main = "Re-start from 10 clusters")
par(opar)

非階層的クラスタリング

サンプル

require(graphics)

# a 2-dimensional example
x <- rbind(matrix(rnorm(100, sd = 0.3), ncol = 2),
           matrix(rnorm(100, mean = 1, sd = 0.3), ncol = 2))
colnames(x) <- c("x", "y")
(cl <- kmeans(x, 2))
plot(x, col = cl$cluster)
points(cl$centers, col = 1:2, pch = 8, cex=2)

## random starts do help here with too many clusters
(cl <- kmeans(x, 5, nstart = 25))
plot(x, col = cl$cluster)
points(cl$centers, col = 1:5, pch = 8)

モデルに基づく手法

  • R
    • パッケージ: mclust
    • 関数: EMclust, mclustBICなど

サンプル

irisBIC <- mclustBIC(iris[,-5])
irisBIC
plot(irisBIC)

subset <- sample(1:nrow(iris), 100)
irisBIC <- mclustBIC(iris[,-5], initialization=list(subset =subset))
irisBIC
plot(irisBIC)

irisBIC1 <- mclustBIC(iris[,-5], G=seq(from=1,to=9,by=2), 
                    modelNames=c("EII", "EEI", "EEE"))
irisBIC1
plot(irisBIC1)
irisBIC2  <- mclustBIC(iris[,-5], G=seq(from=2,to=8,by=2), 
                       modelNames=c("VII", "VVI", "VVV"), x= irisBIC1)
irisBIC2
plot(irisBIC2)

nNoise <- 450
set.seed(0)
poissonNoise <- apply(apply( iris[,-5], 2, range), 2, function(x, n) 
                      runif(n, min = x[1]-.1, max = x[2]+.1), n = nNoise)
set.seed(0)
noiseInit <- sample(c(TRUE,FALSE),size=nrow(iris)+nNoise,replace=TRUE,
                    prob=c(3,1))
irisNdata <- rbind(iris[,-5], poissonNoise)
irisNbic <- mclustBIC(data = irisNdata,
                      initialization = list(noise = noiseInit))
irisNbic
plot(irisNbic)

自己組織化マップ(SOM)

  • R
    • パッケージ: kohonen, som
    • 関数: som, somgrid, plot.kohonenなどなど

サンプル

data(wines)
set.seed(7)

training <- sample(nrow(wines), 120)
Xtraining <- scale(wines[training, ])
Xtest <- scale(wines[-training, ],
               center = attr(Xtraining, "scaled:center"),
               scale = attr(Xtraining, "scaled:scale"))

som.wines <- som(Xtraining, grid = somgrid(5, 5, "hexagonal"))

som.prediction <- predict(som.wines, newdata = Xtest,
          trainX = Xtraining,
          trainY = factor(wine.classes[training]))
table(wine.classes[-training], som.prediction$prediction)

線形回帰分析

  • 回帰分析は最も一般的な手法且つ話題が豊富な手法
  • 線形回帰分析とは
    • 量的データを目的変数とした最も基本的な解析手法(って)
    • 一つの説明変数の場合: 単回帰分析
    • 複数の説明変数の場合: 重回帰分析

サンプル

require(graphics)

## Annette Dobson (1990) "An Introduction to Generalized Linear Models".
## Page 9: Plant Weight Data.
ctl <- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt <- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group <- gl(2,10,20, labels=c("Ctl","Trt"))
weight <- c(ctl, trt)
anova(lm.D9 <- lm(weight ~ group))
summary(lm.D90 <- lm(weight ~ group - 1))# omitting intercept
summary(resid(lm.D9) - resid(lm.D90)) #- residuals almost identical

opar <- par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
plot(lm.D9, las = 1)      # Residuals, Fitted, ...
par(opar)

## model frame :
stopifnot(identical(lm(weight ~ group, method = "model.frame"),
                    model.frame(lm.D9)))

### less simple examples in "See Also" above

非線形回帰分析

  • 線形回帰分析とは
  • R
    • パッケージ: stats, mgcv
    • 関数: nls, glm(一般化線形モデル), gam(平滑化回帰)

線形判別分析

  • 回帰分析との比較
    • 回帰分析: 外的基準が量的データ
    • 判別分析: 外的基準が質的データ
  • フィッシャーの線形判別関数
    • irisの判別などなど
    • 古典的手法
  • 線形判別の注意事項
    • 等分散の制約条件が必要
    • 大量の変数には向かない
  • R
    • 関数:lda
    • predict, cross-validationなど併用する

サンプル

Iris <- data.frame(rbind(iris3[,,1], iris3[,,2], iris3[,,3]),
                   Sp = rep(c("s","c","v"), rep(50,3)))
train <- sample(1:150, 75)
table(Iris$Sp[train])
## your answer may differ
##  c  s  v
## 22 23 30
z <- lda(Sp ~ ., Iris, prior = c(1,1,1)/3, subset = train)
predict(z, Iris[-train, ])$class
##  [1] s s s s s s s s s s s s s s s s s s s s s s s s s s s c c c
## [31] c c c c c c c v c c c c v c c c c c c c c c c c c v v v v v
## [61] v v v v v v v v v v v v v v v
(z1 <- update(z, . ~ . - Petal.W.))

  • 非線形と比べるとあまり使われない

非線形判別分析

  • 非線形判別分析とは
    • 線形判別以外の全部(ぉ)の手法
    • なので、非線形的手法以外にも距離に基づいた判別手法、多数決の判別方法、ベイズ判別方法、機械学習による判別方法も含む。
  • 判別関数による判別分析
    • 二次式に依る判別関数としてRだとqdaがある
  • R
    • パッケージ: MASS
    • 関数: qda

サンプル

tr <- sample(1:50, 25)
train <- rbind(iris3[tr,,1], iris3[tr,,2], iris3[tr,,3])
test <- rbind(iris3[-tr,,1], iris3[-tr,,2], iris3[-tr,,3])
cl <- factor(c(rep("s",25), rep("c",25), rep("v",25)))
z <- qda(train, cl)
predict(z,test)$class
  • 距離による判別分析
    • マハラノビス距離など利用

サンプル

require(graphics)

ma <- cbind(1:6, 1:3)
(S <-  var(ma))
mahalanobis(c(0,0), 1:2, S)

x <- matrix(rnorm(100*3), ncol = 3)
stopifnot(mahalanobis(x, 0, diag(ncol(x))) == rowSums(x*x))
        ##- Here, D^2 = usual squared Euclidean distances

Sx <- cov(x)
D2 <- mahalanobis(x, colMeans(x), Sx)
plot(density(D2, bw=.5),
     main="Squared Mahalanobis distances, n=100, p=3") ; rug(D2)
qqplot(qchisq(ppoints(100), df=3), D2,
       main = expression("Q-Q plot of Mahalanobis" * ~D^2 *
                         " vs. quantiles of" * ~ chi[3]^2))
abline(0, 1, col = 'gray')
  • 多数決による判別分析
    • k-NNなどがある

サンプル

train <- rbind(iris3[1:25,,1], iris3[1:25,,2], iris3[1:25,,3])
test <- rbind(iris3[26:50,,1], iris3[26:50,,2], iris3[26:50,,3])
cl <- factor(c(rep("s",25), rep("c",25), rep("v",25)))
knn(train, test, cl, k = 3, prob=TRUE)
attributes(.Last.value)
  • R
    • パッケージ: e1071, klaRなど
    • 関数: NaiveBayes

サンプル

data(iris)
mN <- NaiveBayes(Species ~ ., data = iris)
plot(mN)

mK <- NaiveBayes(Species ~ ., data = iris, usekernel = TRUE)
plot(mK)

生存分析

  • 生存分析とは
    • イベントが起きるまでの時間とイベントとのあいだの関係の関係に着目した手法
  • 適用範囲
    • 工学: 機械システムや製品の故障
    • 医学分野: 疾患の病気の再発や死亡など
    • 生存分析では故障、破壊、倒産、死亡などのイベントを広義で死亡とみなす
  • R
    • パッケージ: survival
    • 関数: survfitなど

サンプル

leukemia.surv <- survfit(Surv(time, status) ~ x, data = aml) 
plot(leukemia.surv, lty = 2:3) 
legend(100, .9, c("Maintenance", "No Maintenance"), lty = 2:3) 
title("Kaplan-Meier Curves\nfor AML Maintenance Study") 
lsurv2 <- survfit(Surv(time, status) ~ x, aml, type='fleming') 
plot(lsurv2, lty=2:3, fun="cumhaz", 
        xlab="Months", ylab="Cumulative Hazard") 



時系列

  • 目的
    • 時系列データの変動の特徴を捉え、現象の解明と将来の予測、制御をするために利用
  • モデル
    • AR, ARMA, ARIMA, ARFIMA, GARCH, VARなどなど
  • 予測モデルを作る際にはacf, pacf, AIC, スペクトル分析, 場合によっては(複数時系列モデルなど)ccfを事前に確認する必要あり
  • 後は単位根検定
  • R
    • パッケージ: tseries, fracdiff, fseries
    • ar, arma, arima, fracdiff, garch などなど

サンプル

arima(lh, order = c(1,0,0))
arima(lh, order = c(3,0,0))
arima(lh, order = c(1,0,1))

arima(lh, order = c(3,0,0), method = "CSS")

arima(USAccDeaths, order = c(0,1,1), seasonal = list(order=c(0,1,1)))
arima(USAccDeaths, order = c(0,1,1), seasonal = list(order=c(0,1,1)),
      method = "CSS") # drops first 13 observations.
# for a model with as few years as this, we want full ML

arima(LakeHuron, order = c(2,0,0), xreg = time(LakeHuron)-1920)

## presidents contains NAs
## graphs in example(acf) suggest order 1 or 3
require(graphics)
(fit1 <- arima(presidents, c(1, 0, 0)))
tsdiag(fit1)
(fit3 <- arima(presidents, c(3, 0, 0)))  # smaller AIC
tsdiag(fit3)
  • カオス時系列
    • 不規則に変動する時系列データを非線形的に解析する手法
  • R
    • パッケージ: tseriesChaos
    • 関数: embedd

サンプル

library(scatterplot3d)
x <- window(rossler.ts, start=90)
xyz <- embedd(x, m=3, d=8)
scatterplot3d(xyz, type="l")

決定木

  • IF-THENで分岐するtreeを作り描画
    • わりかしよく使われる
  • 分類基準
  • R
    • パッケージ: mvpart
    • 関数: rpart

サンプル

data(car.test.frame)
z.auto <- rpart(Mileage ~ Weight, car.test.frame)
summary(z.auto)
plot(z.auto); text(z.auto)

data(spider)
fit1 <- rpart(data.matrix(spider[,1:12])~water+twigs+reft+herbs+moss+sand,spider,method="mrt")
plot(fit1); text(fit1)
fit2 <- rpart(data.matrix(spider[,1:12])~water+twigs+reft+herbs+moss+sand,spider,method="mrt",dissim="man")
plot(fit2); text(fit2)
fit3 <- rpart(gdist(spider[,1:12],meth="bray",full=TRUE,sq=TRUE)~water+twigs+reft+herbs+moss+sand,spider,method="dist")
plot(fit3); text(fit3)

ニューラルネットワーク

  • 注意事項
    • overfittingしやすい
    • 学習方法によって結果が変わる などなど
  • R
    • パッケージ: nnet

サンプル

# use half the iris data
ir <- rbind(iris3[,,1],iris3[,,2],iris3[,,3])
targets <- class.ind( c(rep("s", 50), rep("c", 50), rep("v", 50)) )
samp <- c(sample(1:50,25), sample(51:100,25), sample(101:150,25))
ir1 <- nnet(ir[samp,], targets[samp,], size = 2, rang = 0.1,
            decay = 5e-4, maxit = 200)
test.cl <- function(true, pred) {
    true <- max.col(true)
    cres <- max.col(pred)
    table(true, cres)
}
test.cl(targets[-samp,], predict(ir1, ir[-samp,]))


# or
ird <- data.frame(rbind(iris3[,,1], iris3[,,2], iris3[,,3]),
        species = factor(c(rep("s",50), rep("c", 50), rep("v", 50))))
ir.nn2 <- nnet(species ~ ., data = ird, subset = samp, size = 2, rang = 0.1,
               decay = 5e-4, maxit = 200)
table(ird$species[-samp], predict(ir.nn2, ird[-samp,], type = "class"))

カーネル法サポートベクターマシン

  • 適用範囲
    • 密度関数の推定、主成分分析、正準相関分析、クラスタ分析、判別分析などなど
  • R
    • パッケージ: kernlab
    • 関数: kpca

サンプル

# another example using the iris
data(iris)
test <- sample(1:150,20)

kpc <- kpca(~.,data=iris[-test,-5],kernel="rbfdot",kpar=list(sigma=0.2),features=2)

#print the principal component vectors
pcv(kpc)

#plot the data projection on the components
plot(rotated(kpc),col=as.integer(iris[-test,5]),xlab="1st Principal Component",ylab="2nd Principal Component")

#embed remaining points 
emb <- predict(kpc,iris[test,-5])
points(emb,col=as.integer(iris[test,5]))
  • R
    • パッケージ: kernlab
    • 関数: ksvm (他にも色々ある)

サンプル

## simple example using the spam data set
data(spam)

## create test and training set
index <- sample(1:dim(spam)[1])
spamtrain <- spam[index[1:floor(2 * dim(spam)[1]/3)], ]
spamtest <- spam[index[((2 * ceiling(dim(spam)[1]/3)) + 1):dim(spam)[1]], ]

## train a support vector machine
filter <- ksvm(type~.,data=spamtrain,kernel="rbfdot",kpar=list(sigma=0.05),C=5,cross=3)
filter

## predict mail type on the test set
mailtype <- predict(filter,spamtest[,-58])

## Check results
table(mailtype,spamtest[,58])


## Another example with the famous iris data
data(iris)

## Create a kernel function using the build in rbfdot function
rbf <- rbfdot(sigma=0.1)
rbf

## train a bound constraint support vector machine
irismodel <- ksvm(Species~.,data=iris,type="C-bsvc",kernel=rbf,C=10,prob.model=TRUE)

irismodel

## get fitted values
fitted(irismodel)

## Test on the training set with probabilities as output
predict(irismodel, iris[,-5], type="probabilities")


## Demo of the plot function
x <- rbind(matrix(rnorm(120),,2),matrix(rnorm(120,mean=3),,2))
y <- matrix(c(rep(1,60),rep(-1,60)))

svp <- ksvm(x,y,type="C-svc")
plot(svp,data=x)


### Use kernelMatrix
K <- as.kernelMatrix(crossprod(t(x)))

svp2 <- ksvm(K, y, type="C-svc")

svp2


#### Use custom kernel 

k <- function(x,y) {(sum(x*y) +1)*exp(-0.001*sum((x-y)^2))}
class(k) <- "kernel"

data(promotergene)

## train svm using custom kernel
gene <- ksvm(Class~.,data=promotergene,kernel=k,C=10,cross=5)

gene


#### Use text with string kernels
data(reuters)
is(reuters)
tsv <- ksvm(reuters,rlabels,kernel="stringdot",kpar=list(length=5),cross=3,C=10)
tsv


## regression
# create data
x <- seq(-20,20,0.1)
y <- sin(x)/x + rnorm(401,sd=0.03)

# train support vector machine
regm <- ksvm(x,y,epsilon=0.01,kpar=list(sigma=16),cross=3)
plot(x,y,type="l")
lines(x,predict(regm,x),col="red")

集団学習

  • 集団学習とは
    • アンサンブル学習とも言われる
    • 決して精度が高いとは言えない分類器の結果から学習を行ない制度が高い分類器の構築を行う
  • 手法
    • バギング、ブースティング、ランダムフォレスト
    • ランダムフォレストは特に新しい。

  • バギング
    • bagging(bootstrap aggregating), 1996年ブライマンによって提案(L.Breiman)
    • 与えられたデータセットをブートストラップというリサンプリング法によって複数の学習データを作成
    • 作成したデータセットに回帰・分析結果を統合、組み合わせることで精度をあげる
    • ブートストラップ、サンプルはそれぞれ独立、学習は並列実行可能.
  • R
    • パッケージ: adabag
    • 関数: bagging

サンプル

## rpart library should be loaded
library(rpart)
data(iris)
names(iris)<-c("LS","AS","LP","AP","Especies")
lirios.bagging <- bagging(Especies~LS +AS +LP+ AP, data=iris, mfinal=10)

## rpart and mlbench libraries should be loaded
library(rpart)
library(mlbench)
data(BreastCancer)
l <- length(BreastCancer[,1])
sub <- sample(1:l,2*l/3)
BC.bagging <- bagging(Class ~.,data=BreastCancer[,-1],mfinal=25, maxdepth=3)
BC.bagging.pred <- predict.bagging(BC.bagging,newdata=BreastCancer[-sub,-1])
BC.bagging.pred[-1]

# Data Vehicle (four classes)
library(rpart)
library(mlbench)
data(Vehicle)
l <- length(Vehicle[,1])
sub <- sample(1:l,2*l/3)
Vehicle.bagging <- bagging(Class ~.,data=Vehicle[sub, ],mfinal=50, maxdepth=5)
Vehicle.bagging.pred <- predict.bagging(Vehicle.bagging,newdata=Vehicle[-sub, ])
Vehicle.bagging.pred[-1]
  • ブースティング
    • 教師付きデータを用いて学習を行ない、学習結果を踏まえ重みの調整を繰り返す
    • 複数の学習結果を求め、結果を統合・組み合わせをすることで精度をあげる
    • AdaBoost(1996)が有名
  • R
    • パッケージ: adabag
    • 関数: adaboost.M1

サンプル

## rpart library should be loaded
library(rpart)
data(iris)
names(iris)<-c("LS","AS","LP","AP","Especies")
iris.adaboost <- adaboost.M1(Especies~LS +AS +LP+ AP, data=iris, boos=TRUE, 
        mfinal=10)

## rpart and mlbench libraries should be loaded
## Comparing the test error of rpart and adaboost.M1
library(rpart)
library(mlbench)
data(BreastCancer)
l <- length(BreastCancer[,1])
sub <- sample(1:l,2*l/3)

BC.rpart <- rpart(Class~.,data=BreastCancer[sub,-1], maxdepth=3)
BC.rpart.pred <- predict(BC.rpart,newdata=BreastCancer[-sub,-1],type="class")
tb <-table(BC.rpart.pred,BreastCancer$Class[-sub])
error.rpart <- 1-(sum(diag(tb))/sum(tb))
tb
error.rpart

BC.adaboost <- adaboost.M1(Class ~.,data=BreastCancer[,-1],mfinal=25, maxdepth=3)
BC.adaboost.pred <- predict.boosting(BC.adaboost,newdata=BreastCancer[-sub,-1])
BC.adaboost.pred[-1]

## Data Vehicle (four classes) 
library(rpart)
library(mlbench)
data(Vehicle)
l <- length(Vehicle[,1])
sub <- sample(1:l,2*l/3)
mfinal <- 25
maxdepth <- 5

Vehicle.rpart <- rpart(Class~.,data=Vehicle[sub,],maxdepth=maxdepth)
Vehicle.rpart.pred <- predict(Vehicle.rpart,newdata=Vehicle[-sub, ],type="class")
tb <- table(Vehicle.rpart.pred,Vehicle$Class[-sub])
error.rpart <- 1-(sum(diag(tb))/sum(tb))
tb
error.rpart

Vehicle.adaboost <- adaboost.M1(Class ~.,data=Vehicle[sub, ],mfinal=mfinal, 
        maxdepth=maxdepth)
Vehicle.adaboost.pred <- predict.boosting(Vehicle.adaboost,newdata=Vehicle[-sub, ])
Vehicle.adaboost.pred[-1]

  • ランダムフォレスト
    • Random Forest(RF)はバギングの提唱者Breimalによって提案された新しい手法
    • 精度、PCの資源節約の面でバギング、ブースティングより優秀
  • R
    • パッケージ: randomForest
    • 関数: randomForest

サンプル

## Classification:
##data(iris)
set.seed(71)
iris.rf <- randomForest(Species ~ ., data=iris, importance=TRUE,
                        proximity=TRUE)
print(iris.rf)
## Look at variable importance:
round(importance(iris.rf), 2)
## Do MDS on 1 - proximity:
iris.mds <- cmdscale(1 - iris.rf$proximity, eig=TRUE)
op <- par(pty="s")
pairs(cbind(iris[,1:4], iris.mds$points), cex=0.6, gap=0,
      col=c("red", "green", "blue")[as.numeric(iris$Species)],
      main="Iris Data: Predictors and MDS of Proximity Based on RandomForest")
par(op)
print(iris.mds$GOF)

## The `unsupervised' case:
set.seed(17)
iris.urf <- randomForest(iris[, -5])
MDSplot(iris.urf, iris$Species)

## Regression:
## data(airquality)
set.seed(131)
ozone.rf <- randomForest(Ozone ~ ., data=airquality, mtry=3,
                         importance=TRUE, na.action=na.omit)
print(ozone.rf)
## Show "importance" of variables: higher value mean more important:
round(importance(ozone.rf), 2)

## "x" can be a matrix instead of a data frame:
set.seed(17)
x <- matrix(runif(5e2), 100)
y <- gl(2, 50)
(myrf <- randomForest(x, y))
(predict(myrf, x))

## "complicated" formula:
(swiss.rf <- randomForest(sqrt(Fertility) ~ . - Catholic + I(Catholic < 50),
                          data=swiss))
(predict(swiss.rf, swiss))
## Test use of 32-level factor as a predictor:
set.seed(1)
x <- data.frame(x1=gl(32, 5), x2=runif(160), y=rnorm(160))
(rf1 <- randomForest(x[-3], x[[3]], ntree=10))

## Grow no more than 4 nodes per tree:
(treesize(randomForest(Species ~ ., data=iris, maxnodes=4, ntree=30)))


アソシエーション分析

  • アソシエーション分析とは
    • POSデータ等から有益な情報を見つける際に活用される
    • 代表的な手法に相関ルール、頻出アイテムなどの手法がある
  • R
    • パッケージ: arules
    • 関数: apriori

サンプル

data("Adult")
## Mine association rules.
rules <- apriori(Adult, 
                 parameter = list(supp = 0.5, conf = 0.9,
                                  target = "rules"))
summary(rules)

サンプル

data("Adult")
## Mine itemsets with minimum support of 0.1.
itemsets <- eclat(Adult,
                  parameter = list(supp = 0.1, maxlen = 15))

2010年07月28日 このエントリーを含むブックマーク

jQueryjQuery(1.4.2)のソースを少し読むなど

  • jQueryのソース読みなど
    • なんとなく実行(ver 1.4.2)
    • 概要掴むくらいの感じで
    • 6000行くらいか(まあ6241行ですけど。。)

for loopのテクとか

  • まあJsでは良く使うテクニックが散見される
    • 例えば下記のinArrayコードのfor ( var i = 0, length = array.length; i < length; i++ )って部分
    • ちなみにjQuery cook bookとかだとeachつかうよりforで回すほうが速いので処理速くする際にはこうする的なことが書いてたり。
    • あと、i++より++iのがブラウザによっては速くなるとか,htmlより[0].innerHTMLのがDOM Updateが速くなるとかいう事もjQuery Cook Bookに書いてたりしたり
inArray: function( elem, array ) {
	if ( array.indexOf ) {
		return array.indexOf( elem );
	}

	for ( var i = 0, length = array.length; i < length; i++ ) {
		if ( array[ i ] === elem ) {
			return i;
		}
	}

	return -1;
},

splitテクとか

  • 何かしらの処理を一気にやるときに使う
    • 一気にevnent bindingしたりするときに使ったり
    • Ajax Event handler attachするときに使ったり
jQuery.each( ("blur focus focusin focusout load resize scroll unload click dblclick " +
	"mousedown mouseup mousemove mouseover mouseout mouseenter mouseleave " +
	"change select submit keydown keypress keyup error").split(" "), function( i, name ) {

	// Handle event binding
	jQuery.fn[ name ] = function( fn ) {
		return fn ? this.bind( name, fn ) : this.trigger( name );
	};

	if ( jQuery.attrFn ) {
		jQuery.attrFn[ name ] = true;
	}
});

ajax関連など(bind bind)

jQuery.each( "ajaxStart ajaxStop ajaxComplete ajaxError ajaxSuccess ajaxSend".split(" "), function( i, o ) {
	jQuery.fn[o] = function( f ) {
		return this.bind(o, f);
	};
});




Ajax周りの関数の実装など

  • 色々な関数があるが、基本的には内部でjQuery.Ajaxを使う構造になっている

global objectとしてのjQueryの登録

  • まあ以下のような感じで
// Expose jQuery to the global object
window.jQuery = window.$ = jQuery;

queueとかdequeueの実装とか

  • overrideで実装(複数あり)

firstとかlastとか

  • selectorとかは単純にeq(0)とかeq(-1)で取得してる
first: function() {
	return this.eq( 0 );
},

last: function() {
	return this.eq( -1 );
},


IE関連の対応実装とか

  • やったらソースコードが長くなってるのはIE周りの対応が原因という節(まあそれほど多いというわけでもなかったり。。)
    • IEの対応部分は基本的にコメントに参考URLが載ってたり
      • memory leak対応が非常に多いのが特徴的(memoryで検索かけると大抵ID対策)
    • まあ、IE以外の部分も結構参考URLコメントに書いとりますがね。

例1: IEでのDOM準備判定とか

// The DOM ready check for Internet Explorer
function doScrollCheck() {
	if ( jQuery.isReady ) {
		return;
	}

	try {
		// If IE is used, use the trick by Diego Perini
		// http://javascript.nwbox.com/IEContentLoaded/
		document.documentElement.doScroll("left");
	} catch( error ) {
		setTimeout( doScrollCheck, 1 );
		return;
	}

	// and execute any waiting functions
	jQuery.ready();
}

例2: Ajaxのアレとか

// Create the request object; Microsoft failed to properly
// implement the XMLHttpRequest in IE7 (can't request local files),
// so we use the ActiveXObject when it is available
// This function can be overriden by calling jQuery.ajaxSetup
xhr: window.XMLHttpRequest && (window.location.protocol !== "file:" || !window.ActiveXObject) ?
	function() {
		return new window.XMLHttpRequest();
	} :
	function() {
		try {
			return new window.ActiveXObject("Microsoft.XMLHTTP");
		} catch(e) {}
	},

IEメモリリーク対策とか

// Prevent memory leaks in IE
// Window isn't included so as not to unbind existing unload events
// More info:
//  - http://isaacschlueter.com/2006/10/msie-memory-leaks/
if ( window.attachEvent && !window.addEventListener ) {
	window.attachEvent("onunload", function() {
		for ( var id in jQuery.cache ) {
			if ( jQuery.cache[ id ].handle ) {
				// Try/Catch is to handle iframes being unloaded, see #4280
				try {
					jQuery.event.remove( jQuery.cache[ id ].handle.elem );
				} catch(e) {}
			}
		}
	});
}

jQuery.fn.extendとかjQuery.extendとか


  • 関数などをする際にはjQuery.fn.extendとかjQuery.extendを利用する
    • 実装は以下のような感じ?(途中省略)
// Give the init function the jQuery prototype for later instantiation
jQuery.fn.init.prototype = jQuery.fn;

jQuery.extend = jQuery.fn.extend = function() {...

謎な部分.ブラウザJs実装に応じてダブってる値の扱いを変更させたりするっぽい?

// Here we check if the JavaScript engine is using some sort of
// optimization where it does not always call our comparision
// function. If that is the case, discard the hasDuplicate value.
//   Thus far that includes Google Chrome.
[0, 0].sort(function(){
	baseHasDuplicate = false;
	return 0;
});

Sizzle.jsの利用

  • CSS Selector enzineのsizzleを内部で利用している

まあざっとソース見ただけだったり.

2010年07月27日 このエントリーを含むブックマーク

Symfony】【PHPSymfony関連資料とか作業工程とかのメモなど

Rails程使うわけじゃないんですけどたまに仕事とかくるんでなんとなくメモメモ。SymfonyDeliciousとかで使われてるらしいですね。他の例はあまり知らないですが。

Symfony関連の資料は普通に刊行されている本がフリーで沢山読めるのが良いですね.一冊各3000円くらい?するので

あと、ヘルパー関数とかはドキュメントとか見るよりは、Symfony本体のソースを見たほうが理解が速いと思いますね。


公式サイト


新しめの資料


少し古めの資料

チートシート(結構新しめ)


まあここらへんの資料読んでおけば大抵なんとかなるかと。



備考: ざっとした作業工程まとめ

まあ個人的な工程もなんとなくメモ用に.

1. テーブル用schemaの更新 + モデル作成

  • vi config/doctrine/schema.yml 変更
  • vi data/fixtures/~ でサンプルデータ追加したり
  • php symfony doctrine:build --all --and-loadとか実行
  • 1テーブル追加するとmodel * 3, filter * 2, form * 2ファイル新規作成されるのでとりあえずgit addなどする
  • vi lib/model/doctrine/~ でモデル関数とか作ったりする

2. routing追加

  • vi apps/frontend/config/routing.ymlとか変更
  • php symfony ccとか実行したり

3. module追加

  • php symfony generate:module frontend nanka
    • 実行すると, action, template, functional testファイルなど作成されたりするのでgit addなどする

4. controller追加

  public function executeNanka(sfWebRequest $request)
  {

  }

5. viewの追加

  • viewはコントローラ名に対応するので対応したファイル名で作成する
  • vi apps/frontend/modules/receive/nanka/template/nankaSuccess.php

後は1~5あたりを繰り返したり、ドキュメントとか本体ソースとか見たりするといった感じかと

2010年07月21日

類似しているプログラミング言語の文法などの比較まとめチートシートが良い感じ

以前見つけた資料。そういやそんなんあったなと久々に検索して探すのに少し手間取ったのでメモ

言語の比較対応で文法覚えられそうなんで便利じゃないかなと


参照: