2019年7月26日金曜日

[Julia] Julia環境下でPythonを利用する

Julia環境下でPythonを利用できるようにする際のメモです.

まずは,ターミナルからJuliaを起動します.
$ Julia
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.1.0 (2019-01-21)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org/ release
|__/                   |

julia> 

"]"を入力してEnterを押して,パッケージの管理モードに入ります(以下のようにカーソルが"pkg>"に切り替わります).

パッケージの管理モードで以下のように入力します.
(v1.1) pkg> add OhMyREPL IJulia PyCall

すると,以下のようにズラズラと表示され,読み込みが始まります.
  Updating registry at `~/.julia/registries/General`

  Updating git-repo `https://github.com/JuliaRegistries/General.git
... 途中省略 ...
  Building PyCall → `~/.julia/packages/PyCall/ttONZ/deps/build.log`
  Building ZMQ ───→ `~/.julia/packages/ZMQ/ItfqT/deps/build.log`
  Building Rmath ─→ `~/.julia/packages/Rmath/BoBag/deps/build.log`
  Building FFTW ──→ `~/.julia/packages/FFTW/loJ3F/deps/build.log`
  Building Arpack → `~/.julia/packages/Arpack/zCmTA/deps/build.log`
  Building GR ────→ `~/.julia/packages/GR/oiZD3/deps/build.log`

終わったら,以下のようにPyCallを使う宣言をします.
julia> using PyCall
[ Info: Precompiling PyCall [438e738f-606a-5dbb-bf0a-cddfbfd45ab0]

julia> 

これで PyCall が使える環境になります(初回は上記のプレコンパイルに多少時間がかかります).

以下は試しにmathライブラリを使用してみた例です.
julia> math = pyimport("math")
PyObject <module 'math' from '/Users/xxx/.julia/conda/3/lib/python3.7/lib-dynload/math.cpython-37m-darwin.so'>

julia> math.sin(90)
0.8939966636005579

julia> math.cos(90)
-0.4480736161291701
なお,Juliaでは標準環境下で sin, cos関数が使えるので,比較のために計算してみると同じ値が得られます.

PyCallの再構築を行うには,以下のように入力します.
julia> import Pkg; Pkg.build("PyCall")
  Building Conda ─→ `~/.julia/packages/Conda/kLXeC/deps/build.log`
  Building PyCall → `~/.julia/packages/PyCall/ttONZ/deps/build.log`

JupyterNotebookで使えるようにするためには,以下のように入力します.
julia> using IJulia
julia> notebook()

これらをJupyter Notebookで実行した結果はGitHubで公開しています.

2019年7月25日木曜日

[R] 正則化パラメータの選択

Lasso推定値は正則化パラメータλの値に依存するため,λの値が異なると当然Lasso推定値も異なってきます.
特に,Lasso推定値の違いは変数選択の結果の違いに直接影響してくるため,正則化パラメータの値の選択は重要です.

前回の投稿(アメリカの犯罪例にLassoを適用)ではアメリカの犯罪データにLassoを適用しました.
その続きとして,用意した正則化パラメータの値ごとに10分割交差検証法を用いて得られた結果を図示すると下図のようになります.
上図において,y軸は誤差2乗和(交差検証法の値)で,x軸は正則化パラメータの対数値になります.赤色の点から伸びている帯は,10分割交差検証法で計算された10個の検証誤差の標準誤差を表します.また,縦の破線は,選択された正則化パラメータの位置を表します.

CVの値が最小になる正則化パラメータの値(λ = 3)に対応するモデルを最適なモデルとして選択しています.

上図を表示するには,前回の投稿(こちら)に引き続いて以下のように入力します.

まずはCVを計算します.
> res.cv <- cv.glmnet(x=X, y=y)

続いて,CV値の推移をプロットします.
> plot(res.cv, xlab="Logarithmic value of regularization parameter", ylab="Regression coefficient")

CV値が最小となる正則化パラメータの値を出力します.
> res.cv$lambda.min
[1] 15.15854

1標準誤差ルールにより選択された正則化パラメータの値を出力します.

> res.cv$lambda.1se
[1] 155.1523


2019年7月24日水曜日

[R] アメリカの犯罪例に Lasso を適用

スパース推定法による統計モデリング (統計学One Point)の 1.3 にある例を実行してみた結果のメモです.

アメリカ合衆国の都市における犯罪データに関して Lasso を適用する.

まずは,Lassoの代表的なパッケージである glmnet を読み込みます.glmnet パッケージのインストール方法はこちら
> ### Lasso
> library(glmnet)
Loading required package: Matrix
Loading required package: foreach
Loaded glmnet 2.0-18

基となるアメリカ合衆国の都市における犯罪データ(crime.txt *)には,以下の6つの項目を50の都市から取得した結果が記載されています.なお,crime.txtはこちらからダウンロードできるデータです.

X1 = total overall reported crime rate per 1 million residents
         人口100万人あたりの犯罪率(目的変数 Y とする)
X2 = reported violent crime rate per 100,000 residents
         人口10万人当たりの暴力的な犯罪率 ← 今回の例では使用しない
X3 = annual police funding in $/resident
         警察の年間資金($ /住民)(説明変数 X1 とする)
X4 = % of people 25 years+ with 4 yrs. of high school
         25歳以上で,高校を卒業した人の割合(説明変数 X2 とする)
X5 = % of 16 to 19 year-olds not in highschool and not highschool graduates.
         16〜19歳のうち,高校に通っていない,もしくは卒業していない人の割合(説明変数 X3 とする)
X6 = % of 18 to 24 year-olds in college
         18〜24歳のうち,大学生の割合(説明変数 X4 とする)

X7 = % of people 25 years+ with at least 4 years of college
         25歳以上で、4年制大学を卒業した人の割合(説明変数 X5 とする)

このデータに対して,X1のデータを目的変数(Y)とし,X3 - X7の5つを説明変数(X1 - X5)として回帰分析を行います.

まずは,データ(crime.txt)を読み込みます.
> crime <- read.table("crime.txt")  # Load crime data

"as.matrix" を用いて行列に変換します.
> crime <- as.matrix(crime)

どのデータを使用するか,また,目的変数,説明変数の設定を行います.
> X <- crime[, 3:7]  # Explanatory variable
> y <- crime[, 1]    # Objective variable

説明変数の標準化,目的変数の中心化を行います.
> X <- scale(X)      # Standardize explanatory variables
> y <- y - mean(y)   # Centralize objective variables
ただし,glmnet パッケージでは,説明変数を標準化しなくても,初期設定で標準化され,係数の推定値は元の尺度に戻して出力されるとのことです.

Lasso推定を行います.
> # Lasso estimation
> res <- glmnet(x=X, y=y)

解パス図の描画を行います.
> # Drawing Solution path
> plot(res, xvar="lambda", label=TRUE, 
+      xlab="Logarithmic value of regularization parameter", 
+      ylab="Regression coefficient", col="black", lwd=2.5) 

正則化パラーメータの値を 20 に固定して,係数の推定値を求めます.
> # Fix value of regularization parameter to 20
> res1 <- glmnet(x=X, y=y, lambda=20)
> res1$beta  # Estimated value of coefficient
5 x 1 sparse Matrix of class "dgCMatrix"
          s0
V3 133.50551
V4 -25.22804
V5  19.45576
V6   .      
V7   .      

上記の例のJupyter Notebookファイルは,GitHubStatical Modeling via Sparse Estimation.ipynbの"Apply Lasso on crime data"で見ることができます.

..-. --- --- - -.--- - . * crime.txt の中身
478 184 40 74 11 31 20
494 213 32 72 11 43 18
643 347 57 70 18 16 16
341 565 31 71 11 25 19
773 327 67 72 9 29 24
603 260 25 68 8 32 15
484 325 34 68 12 24 14
546 102 33 62 13 28 11
424 38 36 69 7 25 12
548 226 31 66 9 58 15
506 137 35 60 13 21 9
819 369 30 81 4 77 36
541 109 44 66 9 37 12
491 809 32 67 11 37 16
514 29 30 65 12 35 11
371 245 16 64 10 42 14
457 118 29 64 12 21 10
437 148 36 62 7 81 27
570 387 30 59 15 31 16
432 98 23 56 15 50 15
619 608 33 46 22 24 8
357 218 35 54 14 27 13
623 254 38 54 20 22 11
547 697 44 45 26 18 8
792 827 28 57 12 23 11
799 693 35 57 9 60 18
439 448 31 61 19 14 12
867 942 39 52 17 31 10
912 1017 27 44 21 24 9
462 216 36 43 18 23 8
859 673 38 48 19 22 10
805 989 46 57 14 25 12
652 630 29 47 19 25 9
776 404 32 50 19 21 9
919 692 39 48 16 32 11
732 1517 44 49 13 31 14
657 879 33 72 13 13 22
1419 631 43 59 14 21 13
989 1375 22 49 9 46 13
821 1139 30 54 13 27 12
1740 3545 86 62 22 18 15
815 706 30 47 17 39 11
760 451 32 45 34 15 10
936 433 43 48 26 23 12
863 601 20 69 23 7 12
783 1024 55 42 23 23 11
715 457 44 49 18 30 12
1504 1441 37 57 15 35 13
1324 1022 82 72 22 15 16

940 1244 66 67 26 18 16

Crime data

Crime data for 50 U.S cities.
The data (X1, X2, X3, X4, X5, X6, X7) are for each city.
X1 = total overall reported crime rate per 1 million residents
X2 = reported violent crime rate per 100,000 residents
X3 = annual police funding in $/resident
X4 = % of people 25 years+ with 4 yrs. of high school
X5 = % of 16 to 19 year-olds not in highschool and not highschool graduates.
X6 = % of 18 to 24 year-olds in college

X7 = % of people 25 years+ with at least 4 years of college


2019年7月23日火曜日

[R] glmnetパッケージをインストール

ターミナルからRを立ち上げます.
$ r

R version 3.6.0 (2019-04-26) -- "Planting of a Tree"
Copyright (C) 2019 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin17.7.0 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

During startup - Warning messages:
1: Setting LC_COLLATE failed, using "C" 
2: Setting LC_TIME failed, using "C" 
3: Setting LC_MESSAGES failed, using "C" 
4: Setting LC_MONETARY failed, using "C" 

glmnetがインストールされてない状態でglmnetを読み込もうとしても,以下のようにそのようなパッケージが存在しないと返ってきます.
> library(glmnet)
Error in library(glmnet) : there is no package called ‘glmnet’

インストールは以下のようにして行います.
> install.packages("glmnet")

すると,以下のようにCRANミラーを選択するよう支持されます.
Installing package into ‘/usr/local/lib/R/3.6/site-library’
(as ‘lib’ is unspecified)
--- Please select a CRAN mirror for use in this session ---
Secure CRAN mirrors 

 1: 0-Cloud [https]                   2: Algeria [https]                
 3: Australia (Canberra) [https]      4: Australia (Melbourne 1) [https]
 5: Australia (Melbourne 2) [https]   6: Australia (Perth) [https]      
 7: Austria [https]                   8: Belgium (Ghent) [https]        
 9: Brazil (PR) [https]              10: Brazil (RJ) [https]            
11: Brazil (SP 1) [https]            12: Brazil (SP 2) [https]          
13: Bulgaria [https]                 14: Chile [https]                  
15: China (Hong Kong) [https]        16: China (Guangzhou) [https]      
17: China (Lanzhou) [https]          18: China (Shanghai) [https]       
19: Colombia (Cali) [https]          20: Czech Republic [https]         
21: Denmark [https]                  22: Ecuador (Cuenca) [https]       
23: Ecuador (Quito) [https]          24: Estonia [https]                
25: France (Lyon 1) [https]          26: France (Lyon 2) [https]        
27: France (Marseille) [https]       28: France (Montpellier) [https]   
29: France (Paris 2) [https]         30: Germany (Erlangen) [https]     
31: Germany (Göttingen) [https]      32: Germany (Münster) [https]      
33: Germany (Regensburg) [https]     34: Greece [https]                 
35: Hungary [https]                  36: Iceland [https]                
37: Indonesia (Jakarta) [https]      38: Ireland [https]                
39: Italy (Padua) [https]            40: Japan (Tokyo) [https]          
41: Japan (Yonezawa) [https]         42: Korea (Busan) [https]          
43: Korea (Gyeongsan-si) [https]     44: Korea (Seoul 1) [https]        
45: Korea (Ulsan) [https]            46: Malaysia [https]               
47: Mexico (Mexico City) [https]     48: Norway [https]                 
49: Philippines [https]              50: Serbia [https]                 
51: Spain (A Coruña) [https]         52: Spain (Madrid) [https]         
53: Sweden [https]                   54: Switzerland [https]            
55: Turkey (Denizli) [https]         56: Turkey (Mersin) [https]        
57: UK (Bristol) [https]             58: UK (London 1) [https]          
59: USA (CA 1) [https]               60: USA (IA) [https]               
61: USA (KS) [https]                 62: USA (MI 1) [https]             
63: USA (MI 2) [https]               64: USA (OR) [https]               
65: USA (TN) [https]                 66: USA (TX 1) [https]             
67: Uruguay [https]                  68: (other mirrors)                

ここでは,40番のJapan(Tokyo)を選択します.
Selection: 40

するとインストール*が始まります.

・・-・ --- --- - -・ --- - ・
* インストール時には以下のように表示されます.
also installing the dependencies ‘iterators’, ‘foreach’

trying URL 'https://cran.ism.ac.jp/src/contrib/iterators_1.0.10.tar.gz'
Content type 'application/x-gzip' length 290575 bytes (283 KB)
==================================================
downloaded 283 KB

trying URL 'https://cran.ism.ac.jp/src/contrib/foreach_1.4.4.tar.gz'
Content type 'application/x-gzip' length 360705 bytes (352 KB)
==================================================
downloaded 352 KB

trying URL 'https://cran.ism.ac.jp/src/contrib/glmnet_2.0-18.tar.gz'
Content type 'application/x-gzip' length 3862714 bytes (3.7 MB)
==================================================
downloaded 3.7 MB

During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
* installing *source* package ‘iterators’ ...
** package ‘iterators’ successfully unpacked and MD5 sums checked
** using staged installation
** R
** inst
** byte-compile and prepare package for lazy loading
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** help
*** installing help indices
** building package indices
** installing vignettes
** testing if installed package can be loaded from temporary location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** testing if installed package can be loaded from final location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** testing if installed package keeps a record of temporary installation path
* DONE (iterators)
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
* installing *source* package ‘foreach’ ...
** package ‘foreach’ successfully unpacked and MD5 sums checked
** using staged installation
** R
** demo
** inst
** byte-compile and prepare package for lazy loading
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** help
*** installing help indices
** building package indices
** installing vignettes
** testing if installed package can be loaded from temporary location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** testing if installed package can be loaded from final location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** testing if installed package keeps a record of temporary installation path
* DONE (foreach)
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
* installing *source* package ‘glmnet’ ...
** package ‘glmnet’ successfully unpacked and MD5 sums checked
** using staged installation
** libs
gfortran  -fPIC  -g -O2  -c glmnet5dpclean.f -o glmnet5dpclean.o
clang -I"/usr/local/Cellar/r/3.6.0_3/lib/R/include" -DNDEBUG   -I/usr/local/opt/gettext/include -I/usr/local/opt/readline/include -I/usr/local/include  -fPIC  -g -O2  -c glmnet_init.c -o glmnet_init.o
clang -dynamiclib -Wl,-headerpad_max_install_names -undefined dynamic_lookup -single_module -multiply_defined suppress -L/usr/local/Cellar/r/3.6.0_3/lib/R/lib -L/usr/local/opt/gettext/lib -L/usr/local/opt/readline/lib -L/usr/local/lib -o glmnet.so glmnet5dpclean.o glmnet_init.o -L/usr/local/opt/gcc/lib/gcc/9/gcc/x86_64-apple-darwin17/9.1.0 -L/usr/local/opt/gcc/lib/gcc/9 -lgfortran -lquadmath -lm -L/usr/local/Cellar/r/3.6.0_3/lib/R/lib -lR -lintl -Wl,-framework -Wl,CoreFoundation
ld: warning: text-based stub file /System/Library/Frameworks//CoreFoundation.framework/CoreFoundation.tbd and library file /System/Library/Frameworks//CoreFoundation.framework/CoreFoundation are out of sync. Falling back to library file for linking.
installing to /usr/local/lib/R/3.6/site-library/00LOCK-glmnet/00new/glmnet/libs
** R
** data
** inst
** byte-compile and prepare package for lazy loading
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** help
*** installing help indices
** building package indices
** installing vignettes
** testing if installed package can be loaded from temporary location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** checking absolute paths in shared objects and dynamic libraries
** testing if installed package can be loaded from final location
During startup - Warning messages:
1: Setting LC_TIME failed, using "C" 
2: Setting LC_MESSAGES failed, using "C" 
3: Setting LC_MONETARY failed, using "C" 
** testing if installed package keeps a record of temporary installation path
* DONE (glmnet)

The downloaded source packages are in
‘/private/var/folders/sd/hzc6sg0j5tj5p340n4h992_r0000gn/T/RtmphxoZOU/downloaded_packages’