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

どらちゃんのポッケ

R・統計・技術メモなど勉強ログ置き場

juliaで線形回帰を試してみて、Rと比較を個人的にしてみる

julia Studio

  • julia Studioでは既にjuliaを抱えている
  • せっかく昨日インストールしたので、MacのHomebrewでインストールしたjuliaを見に行くように設定
  • Tool > external tool > juliaで指定
  • pathは下記な感じ

f:id:sleeping_micchi:20140122012903p:plain

juliaの実行方法

  • hoge.jlファイルがjudiaスクリプト
    • コンソールで対話実行
    • 右上の三角緑で実行
    • include("/Users/judia/hello.jl")

おなじみのHello world

  • println("Hello, World!")

線形回帰をやってみる

元データ

  • http://forio.com/downloads/tutorials/gasoline.csv
  • ガソリンの成分分析
    • 成分1・成分2・成分3の構成比と炭素含有量のテーブル
    • よくわからないが、炭素含有量 = 成分1×a + 成分2×b + 成分3×c + dの回帰式を求める

コード

data = readcsv("gasoline.csv")
x = data[:, 2:4]
y = data[:, 6]
coefs = linreg(x, y)

結果の見方

  • Array{Float64,1}: 102.44 -0.0952646 -0.148308 -0.0831406 が返ってくる
  • 回帰式は下記のようになるらしい
    • 炭素含有量 = 成分1× -0.0952646 + 成分2 × -0.148308 + 成分3 × -0.0831406 + 102.44
  • あれ?統計検定は・・・!?

f:id:sleeping_micchi:20140122013039p:plain

Rとの比較

  • Rで書くとこんな感じ

      calc <- function() {
        startTime<-proc.time()
        d <- read.csv("gasoline.csv",header=F)
        d <- data.frame(d)
        result <- lm(d$V6~d$V2+d$V3+d$V4)
        summary(result)
        return (proc.time()-startTime)
      }   
      calc()
    

回帰式の比較

  • Rでの実行結果

      > lm(d$V6~d$V2+d$V3+d$V4)
    
      Call:
      lm(formula = d$V6 ~ d$V2 + d$V3 + d$V4)
    
      Coefficients:
      (Intercept)         d$V2         d$V3         d$V4  
        102.44000     -0.09526     -0.14831     -0.08314 
    
  • ほぼ一緒ですね。

      > summary(result)
    
      Call:
      lm(formula = d$V6 ~ d$V2 + d$V3 + d$V4)
    
      Residuals:
           Min       1Q   Median       3Q      Max 
      -1.15934 -0.28095 -0.00044  0.36190  1.38077 
    
      Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
      (Intercept) 102.440003   0.682827 150.023  < 2e-16 ***
      d$V2         -0.095265   0.006304 -15.113  < 2e-16 ***
      d$V3         -0.148308   0.038597  -3.843 0.000247 ***
      d$V4         -0.083141   0.012339  -6.738 2.47e-09 ***
      ---
      Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    
      Residual standard error: 0.5331 on 78 degrees of freedom
      Multiple R-squared:  0.8606,    Adjusted R-squared:  0.8552 
      F-statistic: 160.5 on 3 and 78 DF,  p-value: < 2.2e-16
    
  • Rだと、p値とか自由度とかも算出してくれる

実行時間の比較

さて、ここですね。ここ。 ということで、超短いサンプルコードだけどどうかなーってはかってみました。

  • R

      > calc()
         ユーザ   システム       経過  
       0.004      0.000      0.004 
    
  • julia(IDEから実行1回目)

      julia> LinearRegression
      1.5374410152435303
    

ん??

  • julia(IDEから実行2回目)

      julia> LinearRegression
      0.0019729137420654297
    

ほう

f:id:sleeping_micchi:20140122013236p:plain

  • julia(コンソールから実行)

      michi$ julia LinearRegression.jl
      2.358590841293335
    

あれ??

  • 結果
    • juliaの2回目 > R > julia(IDEから実行1回目) > julia(コンソールから実行)となりました
    • 2回目の方が本当の実力なら確かに早いな
    • 1回目が遅いのは、コンパイル時間のせいっぽい

本当に短いサンプルを書いた個人的感想(確かめないといけない部分含む)

  • juliaの統計解析パッケージでp値などが算出されないのか?
    • もし算出されないとなると統計解析はRの方がよいと思う
      • その辺の値は、自分で算出しろよって話もある・・・
      • でないはずはないと思うんだけどな。
    • 結果の出力もラベルなしの配列で返ってくる
      • アドホックな分析で使う場合、ラベルがついていないとつらい。
  • 個人的にRに慣れすぎてしまっているから、julia使いにくいなーって思う・・・
  • IDE
    • RStadioの完成度がかなり成熟してきている感があり、juliaStudioの方はまだきついなって感じ
    • Ijuliaの方は試していないが、RStadioの方が簡単に扱えると思っている
  • ライブラリ
    • Rのほうが長いってのもあって、ライブラリは充実している。関数にぶち込むだけなら大抵の分析は知識なくてもできてしままう
    • juliaはRのパッケージにくらべると充実度は低いけど、AWSとかDokerとかあって、エンジニア臭がRよりも高いって第一印象
      • juliaにもMCMCとか、GLMとかのパッケージはあるっぽいから基本的な統計処理を行う分には問題ないと思う
  • 早いっていうけど、早さを実感できなかった
    • まだちゃんとドキュメントを読んでない(すいません)なので、どこが早いとか初回起動が遅いとか、得意な処理・苦手な処理とか、もろもろのことを押さえていないが・・・
  • 学習コスト、ライブラリ、環境などなどを考えた場合、Rの方がまだしばらく存在感を発揮しそうな感じ
    • ドキュメント・パッケージ・認知度とかがあがって技術が枯れてきたらRの置き換えが一気に進みそうな気もする
      • Rと似た感じのコードの書き方なので、乗り換えは容易か?
    • そもそもRを使うような環境でパフォーマンスを真剣に求めないってのもあり
  • juliaを触るなら、Python触りたいのでjuliaはしばらく様子見かなー?