読者です 読者をやめる 読者になる 読者になる

horioの雑記帳

データ分析屋さんが、気の向いた事を色々メモってる雑記帳。

虹色の雲でRを動かす

※この記事は R Advent Calendar 2015 - Qiita の 12日目の記事です。

最近、何故か色々なクラウドで分析のお仕事をしてます(泣)。私の頭は、それぞれの雲の色が入り乱れて真っ黒に?なりそうなので、備忘録がてらに各クラウドでのRの使い方をまとめてみました。

以下の、基本方針や注意事項や

  • 便利なモノは使うに限る。昔はイチから作る気概があった気もしましたが、最近はすっかり易きに…。
  • 各サービスでのアカウントの作り方、料金体系、無料枠の範囲等については、各サービスでの最新情報をご参照ください。
  • 色の表記は私の適当です。公式ではこう言っている、なツッコミはご教示頂ければ幸いです。
  • この色は虹にないぞ、なツッコミは無視しますw
  • (言い訳ですが)12日目を確保した数日後、ドンピシャで週末を全部潰す地方遠征が決定してしまいました…。後で、加筆UPDATEします。

AWS(オレンジ)

  • 後で、19日目の記事にリンクします。

Azure(空色)

最近、見た目や操作方法が大幅に変わり、私は只今大混乱中です…。
なお以下Azureの記述の一部は、既にぞうさん通信(2015/12/7)でも挙げられております。

Rの使い方その1:AzureのVMを使う

極力楽をするとの方針に従えば、Microsoft's new Data Science Virtual Machineで紹介されているとおり、Microsoft Data Science Virtual Machineを使うに限ります。具体的な手順は、以下記事をご参照下さい。
azure.microsoft.com

なお、用意されているものは以下のとおりです。普通に仕事で分析をするには、ほぼ困らなさそうな品揃えです。

Rの使い方その2:Azure Machine Learning(ML)のRノードを使う方法

手軽に高性能のアルゴリズムを使いたい*1、自作APIを作りたい、のときの有力な選択肢と思います。
(少々脱線ですが)12/10に、Azure Stream AnalyticsからAzure MLを操作出来るようになったとの事です。なんかにかこつけて、遊んでみたいと思います。
Azure ML Now Available as a Function in Azure Stream Analytics - Machine Learning - Site Home - TechNet Blogs

ここでは、MLでのRノードでRコードを記述する際の注意点を挙げます。
(後で色々加筆します。長くなったならば、ページを分割します。)

  1. INPUT/OUTPUTは、data.frame形式のみしか受け付けない
  2. Azure ML上で用意されていないpackageを使う際は、依存pakeageも含めて全て、事前にzip形式でUPする必要がある
    • 用意されているpakeage一覧を取得する方法を書く
    • install.packageでのpathの指定方法がややトリッキー(なので、何かしたら上手く成功したのですが、肝心の「何をしたのか」のメモを紛失してしまい、只今大弱り。)
  3. この記述によると、並列実行が出来ない、との事。Azure MLでは、折角A8インスタンスIntel Xeon E5-2670 8コア@2.6 GHz + DDR3-1600 MHz 56 GB)を使っているのに。

参考

Google Cloud Platform(水色+赤+黄色)

GCPでRを立てた経験はないです。易きに流れる、の基本方針に沿うと、こんな記事が良さげに見えました(が、当然実際の検証はしておりません。)
markedmondson.me

IBM BlueMix(水色+黒)

(後で色々加筆します。長くなりそうなので、恐らくページを分割します。)

BlueMixでのPaaSの構成方法

現時点(2015年12月)では、以下の2(+1)通りあります。

  1. ボイラーテンプレートの利用
  2. ランタイム or コンテナ に、APIを乗せる*2
    1. 自作のランタイム or コンテナをUPする
R以前にハマったドツボ

Rの使い方その1:dashDBを使う

  • デフォルトで、dashDB(DB2 ?)+ R Studioが用意済。起動はマウスクリックのみで可能
  • DBとのやりとりは、ibmdbRパッケージを使う

Rの使い方その2:node-REDから使う

  • IoTのセンサーデータ解析なら、これをやりたい!
  • 基本戦略は、dashDBのインスタンスに、POSTコマンドを打つ。結果はDBに保存する
  • DBとのやりとりは、ibmdbRパッケージを使う
注意点
  • sudo権限が無いので、依存関係でsudoを要求するパッケージ(例えばdevtools)は利用できない。痛い
  • 不揮発データは、DBに格納するしかない。node-REDからgetwd()すると、毎回pathが変わるという謎
  • DBのデータが1Gを越えると有料枠に掛かる。dashDBは2015年12月現在で月5250円~、と結構お高い

*1:回帰問題で、RでAzure MLに勝負を挑んだこともありますが、結局勝てませんでしたorz。

*2:ボイラーテンプレートは、既にAPIを乗せているランタイム or コンテナ

ggplot2で一度に複数ページのグラフを出力したい

データ分析の基本は可視化!で、1000変数のオーダーでもプロットして確認したくなりません?な。
もちろん、1000変数のオーダーだと、1ページにまとめられる訳ないので、1ページ1プロットをPDF出力して、パラパラ漫画式のチェックで*1

過日某氏より、上記な事をしたいけど、ggplot2だと複数ページの出力が出来ないと言われたのでメモ。ポイントは明示的にprintをすること。

library(ggplot2)
data(iris) # サンプルデータ

pdf("hoge.pdf")
    for (i in seq_along(names(iris))){
        p <- ggplot(iris, aes(x=Sepal.Length, y=iris[,i])) +
             geom_point() + ylab(names(iris)[i])
        print(p) # ココがミソ
    }
dev.off()

残された疑問

上記、Rっぽくないコードな。applyなりをかます方法は無いだろうか?と。

*1:膨大になると確認の見落としが当然出るので、平均・分散、等々の適当な数値を算出し、変数のスクリーニングやチェックの優先順位つけも併用します。

DBで重複行を削除する

たまに突然必要になるが、忘れそうなのでメモ。OLAP関数(Window関数、分析関数とも)一発で済むのがご利益。確認してないけど、OLAP関数が使えるDBならどれでも行けるハズ。

重複行が発生すること自体、絶対DB設計が間違ってるよね、とのヤジには、ワタシは心中では大いに賛同するものの、口にする勇気はないです…。

重複行削除のSQL

-- テーブルの用意
create table OLAPDeleteT(Val1 integer,Val2 integer);
insert into OLAPDeleteT values(1,1),(1,2),(1,2),(1,2),(1,2),(2,1),(2,1);

-- 以下が本体
delete from (select Row_Number() over(partition by Val1,Val2) as rn
               from OLAPDeleteT)
where rn > 1;

上記SQLのコピペ元。なおこのリンク先は、SQL中級者以上なら、一度はじっくり読みこむ価値がある、オススメページ。
DB2の分析関数の使用例

使う際の留意点

  • ほぼ間違いなくテーブルスキャンになるので遅い*1。レコード数が多い場合は*2、deleteのためだけに、一旦仮に全列にindexを張り、軽くででも統計情報を取り直す事も、十分に合理的な場合もある。
  • (上と関連するが)DBのソートヒープの大きさの確認。行数・列数が多いと、テーブルスキャンなのでそれなりにヒープを食うので、他人の作業の足を大いに引っ張るとか、クエリが実行できない、等々が起こり得る。

*1:全列index張っているなら別だけど、まず有り得ないでしょう…。

*2:個人的には、100万件以上なら警戒する。

半角カタカナを全角カタカナに置換する

とのやりとりを元に、少々改良してみた次第。
なお私は、半角カタカナでもExcelの束でも何でも、やれと言われればやらざるを得ないので、半角カタカナを扱う職に…といった躊躇は何もないとか。

なおNippon::zen2hanのヘルプにはこんな記述が。当然、多いに賛同するのですが、私には向き合わねばいけない現実が、、、とか。
http://cran.r-project.org/web/packages/Nippon/Nippon.pdf

The targeted zenkaku characters are shown with zenkaku constant built into 
Nippon package: only alphabets and numbers. Katakana is not the target of 
zen2han because the halfwidth Katakana is rather a troublemaker.

ワイン方程式の論文とデータを見てみた

タイトルの通りです。なぜ調べたのか?は自分でも疑問。

ワイン方程式とは?

その数学が戦略を決める (文春文庫)

その数学が戦略を決める (文春文庫)

の冒頭で紹介されている、

ワインの品質 = 12.145 + 0.00117×冬の降雨量
             + 0.0614 ×育成期平均気温
             - 0.00386×収穫期降雨量

の式を指すことにします。
他の解説としては、これを参照。

言いたき事

元の論文とデータを探してRで再現したら、上式が間違っていることが分かった。
正確に書くと、こんな式に。

ln(あるビンテージの平均価格/1961年ビンテージの平均価格)
= -12.145 + 0.00117×冬の降雨量
          + 0.614 ×育成期平均気温
          - 0.00386×収穫期降雨量
          + 0.0239 ×1983年基準でのワインの熟成年数

上記の式と比べると、切片項の符号が違うし、ワインの熟成年数の項が抜けてますね。。。
【1/20追記】育成期平均気温の桁も違ってますね。私の転記ミスを修正しました。

現論文&原データのありかと、式の計算証跡は以下の通り。

間違いの発生原因はどこか?

既にまとめられている方がいました。このポストを書きながら、初めて気づいた次第な…。

つまり

  • 原著の英語版で、既に引用が間違っていて、大々的に拡散
  • 日本語訳にも、そのまま間違いが引き継がれた
  • 間違った式が、日本でも各所に拡散した

うん。引用する際には原典に当たるって大事、って再確認。

2015/1/17追記

上記ココログさんより前に、@berobero11 さんが書かれていました。私の調査不足ですみませんでした。ご指摘、ありがとうございました!

beroberoさんのブログポスト:
アッシェンフェルターのワイン価格予測式をトレースしてみた

DB2 V10.5の権限

 忘れそう&案外資料がないのでメモっておく。地味なネタで誰得なのは承知済。

話の前提

ORACLEを念頭にした際の留意点

本題:V10.5の権限関係

 DB2の権限は、V9.7以降大幅に変わっており、管理者が3人必要。ポイントは、それぞれの出来る事が、互いに背反であること
例えば、SYSADMはDBを作成できるものの、作成したDBの中身に対してはLOAD(他DBで言うBULK INSERT)が打てないとか。

  • SYSADM:インスタンスレベルでの操作に関する管理者。例えば、DBの作成など
  • SECADM:与えられたDBに対する、セキュリティ関連の管理者。例えば、GRANT, REVOKEなど
  • DBADM:与えられたDBに対する、DBレベルでの操作関連の管理者。例えば、DROP/CREATE TABLE, REORG, RUNSTATS, LOAD、等々。

 と文字で書くより、以下の図を見た方が早い。
f:id:horihorio:20150109233329p:plain
f:id:horihorio:20150109233643p:plain

 なおこの図は、以下DB2 V10.5公式マニュアル(日本語訳)の、PDFでのページ数p.42とp.44よりコピペ。
http://public.dhe.ibm.com/ps/products/db2/info/vr105/pdf/ja_JP/DB2Security-db2secj1051.pdf

留意点

 DBを作る際にIDが3つ必要。同一IDで、異なる系統の権限を兼任しようとしたら怒られた。そのため、DBADMで投げたクエリの実行状況を取得するために、そのままSYSIBMADMスキーマが触れずに大騒ぎが発生とか。業務では、運用設計や書類、各所への調整などで結構高いハードルかも。

(ついで)~V9.5までと、V9.7との権限の比較

http://www.ibm.com/developerworks/jp/data/products/db2/pdf/db2v97_security.pdf
これの、PDFでのページ数pp.5-6を参照。

取得データからデータマートを作ろうとしてみた

※この記事は R Advent Calendar 2014 - Qiita の 19日目の記事です。

 前編RでXBRLデータを取得してみた - horioの雑記帳では、XBRL形式での有価証券報告書を、zip形式でローカルに保存しました。
 後編の本ポストでは、保存したzipから、表形式での財務諸表データを作成してゆきます。

前・後編全体を通した手順

 前編の再掲ですが、このような流れです。
f:id:horihorio:20141215223656p:plain

(前編)RでXBRLデータを取得してみた → ここRでXBRLデータを取得してみた - horioの雑記帳
(後編)取得データからデータマートを作ってみた → 本ポスト

プログラム本体

 今回も、著者より献本頂いた本のソースが下敷きになっております。下記リンクより、書籍に掲載されいるサンプルコードのDLも可能です。

1. 書籍のサンプルコードからの改変点

  • サンプルコードでの「サンプルインスタンス」を、解凍したzip下の"*.xbrl"に置き換える
  • 日本語名称は、横持ちさせてから最後に紐付ければ良いので省略
  • 連結・単体の両方の要素を取得するように変更。理由は以下の通り
    • 単体=連結の企業の場合、連結の要素が全部NULLの場合があるから
    • 両方でも一方でも計算負荷はほぼ一緒、かつ後で削除するのは簡単なので
  • 処理のボトルネックはCPUなので、福島氏の著書や大仏様のお告げを参考に、並列処理に変更



2. 使い方

以下条件を満たせば、コピペでそのまま動きます。なお今度は、外部サービスへの迷惑はかからないので一安心です。

  • 前編でDLしたzipが、プログラムを置いたディレクトリ直下の"data"フォルダー下に全てあること

3. プログラム

INPUT
  • 前編でDLしたzip
OUTPUT
  • "definfo.txt" -> ID、名称、連結・単体の区分、会計期間。中身はこんな感じ
"ED2013062600002"	"株式会社トランスジェニック"	NA	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600002"	"株式会社トランスジェニック"	NA	"CurrentYearConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600003"	"コムテック株式会社"	NA	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600003"	"コムテック株式会社"	NA	"CurrentYearConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600004"	"株式会社DTS"	NA	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600004"	"株式会社DTS"	NA	"CurrentYearConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600005"	"株式会社フライングガーデン"	NA	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600006"	"㈱寺岡製作所"	"TERAOKA SEISAKUSHO CO., LTD"	"CurrentYearConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600006"	"㈱寺岡製作所"	"TERAOKA SEISAKUSHO CO., LTD"	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600007"	"(株)大戸屋ホールディングス"	"OOTOYA Holdings Co., Ltd."	"CurrentYearConsolidatedDuration"	"2012-04-01"	"2013-03-31"
"ED2013062600007"	"(株)大戸屋ホールディングス"	"OOTOYA Holdings Co., Ltd."	"CurrentYearNonConsolidatedDuration"	"2012-04-01"	"2013-03-31"
  • "FI.CurrentYear.txt" -> DLした縦持ち財務諸表。中身はこんな感じ
"ED2013062600002"	"jpfr-t-cte_IssuanceOfNewSharesExerciseOfSubscriptionRightsToSharesNA"	"CurrentYearNonConsolidatedDuration"	"JPY"	"2186000"
"ED2013062600002"	"jpfr-t-cte_NetIncomeNA"	"CurrentYearNonConsolidatedDuration"	"JPY"	"18877000"
"ED2013062600002"	"jpfr-t-cte_NetChangesOfItemsOtherThanShareholdersEquityNA"	"CurrentYearNonConsolidatedDuration"	"JPY"	"-4391000"
"ED2013062600002"	"jpfr-t-cte_TotalChangesOfItemsDuringThePeriodNA"	"CurrentYearNonConsolidatedDuration"	"JPY"	"16672000"
"ED2013062600002"	"jpfr-t-cte_IssuanceOfNewSharesExerciseOfSubscriptionRightsToSharesCAP"	"CurrentYearNonConsolidatedDuration"	"JPY"	"1093000"
"ED2013062600002"	"jpfr-t-cte_TotalChangesOfItemsDuringThePeriodCAP"	"CurrentYearNonConsolidatedDuration"	"JPY"	"1093000"
ソースコード:"02.perseXBRL.R"
rm(list=ls())
invisible(gc()); invisible(gc())

##### library #####
library("XBRL")
library("foreach")
library("doParallel")

##### prerequirement #####
setwd("./data/")
# unzip対象リスト取得
str.filelist <- dir(pattern = "*zip")

# 並列処理準備
n.cores <- detectCores()
cl <- makeCluster(n.cores, type = "PSOCK")
registerDoParallel(cl)

##### main #####
rm(f.export.persed)
f.export.persed <- function(str.filelist, i){
### zip2xbrl
  # zip内ファイル名取得
  wk.file.list <- unzip(zipfile = str.filelist[i] ,list = TRUE)
  # zipファイル解凍
  unzip(zipfile = str.filelist[i])
  unlink("XbrlDlInfo.csv") # 不要なので即消去
  
  # 拡張子が"xbrl"のファイル名を取得し整形
  wk.file.xbrl <- wk.file.list[grep("xbrl", wk.file.list[,1]), 1]
  wk.file.xbrl <- paste0("./", wk.file.xbrl)
  ### xbrlをparse ### この一行が非常に重い! ###
  xbrl.parsed <- xbrlDoAll(wk.file.xbrl) 
  
  # cleanup
  # 解凍フォルダを削除
  wk.dirname <- unlist(strsplit(wk.file.xbrl,"/"))[2]
  unlink(wk.dirname, recursive = TRUE)
  # 一時変数を全削除
  rm(list = ls(pattern = "wk*"))

### 文書情報、財務諸表情報を抽出
  dat.facts <- xbrl.parsed$fact
  # 内容(=facts) の文字コードを,UTF-8 からShift-JIS に修正
  dat.facts$fact <- iconv(x=dat.facts$fact,from="UTF-8",to="Shift-JIS")
  
### 文章情報
  # 開示者名称を抽出(日英両方)
  target <- grepl(pattern="^DocumentInfo",x=dat.facts$contextId)
  dat.Presenter.Info <- dat.facts[target,]
  dat.Presenter.Info <- dat.Presenter.Info[grepl(pattern="^jpfr-di_EntityName", x=dat.Presenter.Info$elementId), "fact"]
  
  # 期首・期末日を抽出(連結・単体)
  dat.contexts <- xbrl.parsed$context
  target <- grepl(pattern="^CurrentYear",x=dat.contexts$contextId)
  dat.contexts <- na.omit(dat.contexts[target, c("contextId","startDate","endDate")])
  
  # 解凍元ファイル名をKeyとして付与
  dat.definfo <- cbind(sub(".zip","",str.filelist[i]), dat.Presenter.Info[1], dat.Presenter.Info[2], dat.contexts)
  # ファイルに追記出力
  write.table(dat.definfo, "../definfo.txt"
              , sep = "\t", col.name = FALSE, row.names = FALSE, append = TRUE)

### 財務諸表データ
  # 当期要素を抽出(連結・単体)
  target <- grepl(pattern="^CurrentYear",x=dat.facts$contextId)
  dat.FI.CurrentYear <- dat.facts[target,]
  # 解凍元ファイル名をKeyとして付与
  dat.FI.CurrentYear <- cbind(sub(".zip","",str.filelist[i]), dat.FI.CurrentYear[,c(1:4)])
  # ファイルに追記出力
  write.table(dat.FI.CurrentYear, "../FI.CurrentYear.txt"
              , sep = "\t", row.names = FALSE, col.name = FALSE, append = TRUE)
  
### cleanup(これをしないとメモリーが逼迫)  
  rm(list = ls(pattern = "^dat|target|^xbrl"))
  invisible(gc()); invisible(gc())

### 経過出力
  print( c(i, format(Sys.time(), "%m/%d %T")) )
} 

##### main #####
foreach(i = 1:length(str.filelist)
        , .export = ls(envir=parent.frame())
        , .packages = c("XBRL") ) %dopar% {f.export.persed(str.filelist, i)}

stopCluster(cl)

上記プログラムの問題点

 以下4点、個人的な対応コスト順に挙げます。一部の問題は、XBRL仕様書を読めば解決しそうな気がします*1

1. 今年以降に提出される、次世代EDINETタクソノミに未対応

 次世代EDINETタクソノミは何か?は、金融庁HPをご参照下さい。
 重要な点は、上記プログラムは過去分のみにしか通用しないことです。ですので、次世代EDINETタクソノミ用のプログラムを用意する必要があります。なんか悪い予感はしており…。
 また、今年の提出物を眺めると、過去形式で提出された模様なのが散見されます。なので、zipを解凍した結果を見てから、適用するのは 上記ソース/新規のプログラム を判定するロジックも考えないといけません。ただ、数ファイルを実際に眺める限りだと、ロジックと実装は簡単そうな気がしております。

2. IFRS国際会計基準)未対応

 日本でのIFRS適用第1号:日本電波工業(6779)の2013年については、上記ソースではダメでした。1.と同じ話なのかも、も含め、まだ調べておりません。

3. 会計基準が分からない

 現在日本では、日本/米国/IFRSの3基準から選択可能ですが、XBRLのみの情報で、会計基準を判断する情報がまだ分かっておりません。
 なお、どの会計基準なのか?の情報は非常に重要です。例えば売上高については*2、以下リンクの住友商事JTを見ると、会計基準だけで数倍のオーダーで数字が変わります*3。よって、考慮しない訳にはいきません。

4. xbrl::xbrlDoAllが激重

 上記プログラムの処理時間は、9割程度がxbrl::xbrlDoAllになります。CPU勝負なので、AWSのc3.largeで2並列処理をかけてみると、1ファイル平均で20秒程度になります*4東証1/2部、JASDAQ、マザーズで上場企業は約3400社あるので、全企業1年分ならば3/4日ですか。
解決策としては、次の2点を考えてます。

  • AWSにお金を積む

 高級インスタンス借りれば?と思い実験もしましたが、vCPU/ECUが上がるほど処理は早くなるものの、やはり時間短縮は収穫逓減でした。よってお財布との天秤の結果、コンピューティング最適化シリーズで最安のc3.largeにしました。
ついでに無料枠のt2.microも実験しましたが、1日の処理数が100との計算でした…。精神衛生上、やらないことにしました。

  • xbrlDoAllを使わず、XMLで実装する

 取得したい項目が決まっているのだから、この筋でも行けるかも?とも思っています。次世代EDINETタクソノミ対応の際に、気が向けば挑戦するかもしれません。

と、12/17夜時点で下書きしてましたが

 12/18朝の通勤電車内で衝撃が。

 arelleは知らなかったですし、Elasticsearch面白そう、金融データで遊べるのは面白そう、とかとか思ったのと同時に、心からポキッと音がしたのでした。。。
 ただ、よくよく記事を読むと、情報の取得源は有料なのでしょう(多分)。XBRLでどんな情報が不足とか変とか、出来る事を何とかやってみよう、と自分に言い聞かせるのでした。

5. FI.CurrentYear.txt(縦持ち財務諸表)の横持ち変換ソース、各種結合を書いていない

 全上場企業約3400社を1年間ならば、ファイルサイズは120M程度の計算です。この程度ならばRでも出来るかも?と勝手に思っています。出来ないならば、Perl/Pythonの召喚とか。
 各種の出力の結合については、id([id名].zip が随所に登場)が各出力データにあるので、その項目を用いて各種結合を行う予定です。

感想

  • Elasticsearchの記事で、心がぽっきり折れかかっております
  • こういうときは、オレ何やりたかったのだっけ?を今一度考え直すことでしょうか。年末、旅に出ます。。。

Disclaimer

本記事は個人的な勉強のために作成したものです。本情報の内容については万全を期しておりますが、その内容を保証するものではありません。これらの情報によって生じたいかなる損害についても、筆者は一切の責任を負いません。
なお、本資料における意見、見解等はすべて筆者の個人的なものであり、筆者の属する組織の意見、見解等ではないことをご了承ください。

*1:これまで読まずにやってきたのかい!とのツッコミは正論です。はい。

*2:適当にググった結果なので数字が古いです

*3:数倍も違う理由は、リンク先のとおりです。商品原価の大半が税金になる商材の場合には良くある話です。

*4:処理時間の計算は、私の目算w