でたぁっ 感動と失敗の備忘録

データ解析を担当することになったサラリーマンの備忘録

PostgreSQLのcopyではまった

PostgreSQL SQL

PostgreSQLのcopyではまったのでメモしておく。

環境はWindows7 64bit/32bit、PostgreSQLは9.4.4。下記のような商品マスターっぽいtsvデータをcopyすると

id_item name_item
111111 hogehoge ¥200
222222 fugafuga ¥300 

copy table_item from 'item.tsv' encoding 'sjis';


ERROR: invalid byte sequence for encoding "UTF8": 0x80
SQLステート:22021
コンテキスト:COPY table_item, line 1: "111111	hogehoge \200"

ということは、¥ → ¥¥ にすれば良いんだね。

id_item name_item
111111 hogehoge ¥¥200
222222 fugafuga ¥¥300

みたいにすれば当然OK

しかし、今更データいじりたくないので調べるとこんな記事がstackoverflow.com

with で csv形式のコピーにしといて、デリミタでタブを指定です。

copy table_item from 'item.tsv' with csv delimiter e'\t' encoding 'sjis';


クエリーは、成功しました:  2 行の影響があり, 実行時間は、11 ミリ秒でした。

これでうまくいくのですが、まだ問題が！
データ中に引用符「"」があると。。。

id_item name_item
111111 hogehoge" ¥200
222222 fugafuga" ¥300 

Query returned successfully: one row affected, 15 ms execution time.

んっ？？？　エラーにはならないけど。。。
結果1行しかinsertされてない。しかもitem_nameに2行目が入っちゃってます。

postgres=# select * from table_item;
 id_item |       name_item
---------+-----------------------
 111111  | hogehoge \200\r      +
         | 222222  fugafuga \300
(1 行)

データによってはこんなエラーが出るケースも

ERROR:  unterminated CSV quoted field

引用符が悪さしているので、データの中に発生しそうにないバックスペース(x08)を引用符として指定してみた。

copy table_item from 'item.tsv' with csv delimiter e'\t' quote e'\x08' encoding 'sjis';

とりあえず、これでタブ区切りのデータcopyでエラーが出ることはなくなった。
こんな対応で良かったのだろうか。もっとエレガントな解決方法があるかもしれない。

PostgreSQLを使ってみたい、SQLを学んでみたい方へ

　私は1995年頃から業務でSQL書き続けているが、この書籍は入門者はもちろん、ある程度知っている方にもお薦め。最近はWindow関数が充実して便利になったものだ。

10年戦えるデータ分析入門 SQLを武器にデータ活用時代を生き抜く (Informatics &IDEA)

10年戦えるデータ分析入門 SQLを武器にデータ活用時代を生き抜く (Informatics &IDEA)

networkD3

R ネットワーク分析

　RStudioから出たhtmlwidgets for Rって可視化のデパートっぽいな。

他にnetworkD3(ネットワーク図)とdatatables(表)がある。興味があるネットワーク図を試してみた。

networkD3とは

　networkD3は、RからD3でネットワーク系グラフを簡単に作成するためのパッケージ。以前ブログにも書いたd3Networkがベースになっているっぽい。(パッケージ名まぎらわしく間違えそうになる)

# シンプルなネットワーク図を描いてみる
install.packages("networkD3")
library(networkD3)
library(igraph)
simpleNetwork(get.data.frame(barabasi.game(200)))

他にも、forceNetwork、sankeyNetwork、treeNetworkを描くことが出来る。動くものはRPubsにUPした。

他に、HTMLファイル出力

library(magrittr)
simpleNetwork(networkData) %>%
  saveNetwork(file = 'Net1.html')

出来たり、Shinyでも動く。

shiny::runGitHub('christophergandrud/networkD3-shiny-example')

サンプルソースはnetworkD3/inst/examples/shiny at master · christophergandrud/networkD3 · GitHubにあるので、ローカルに保存して

library(shiny)
runApp(appDir="保存先フォルダ")

Shiny - Single-file Shiny apps用のコードはchristophergandrud/networkD3-shiny-example · GitHubにサンプルapp.Rがある。

library(shiny)
runApp(appDir="app.R保存先フォルダ")

ui.Rとserver.R、app.R が一緒に保存されている場合は、ui.Rとserver.Rの方が優先されるようだった。

cyREST API で R と Cytoscape の連携を試す

R Cytoscape API ネットワーク分析

　cyREST APIを使うと、RとCytoscapeを連携させることが出来る。
Cytoscapeは機能を拡張することが出来るアプリ(プラグイン)が数多く提供されており、cyRESTもその中の1つ。

1.CytoscapeにcyRESTをインストールする

　Cytoscapeを起動、Apps→App ManagerでApp Managerを開く。アプリの一覧の中からcyRESTを選択し右下にあるInstallボタンをクリック。インストール済アプリは、Currently Installedタブで確認出来る。
f:id:deta:20141227152704j:plain
　正常に動作していれば、ブラウザに

http://localhost:1234/v1

を入力すると、結果(PCのコア数、メモリサイズなどのJSON)が表示されるはず。他にもJARファイルをダウンロードしてインストールなど、いくつかの方法かある。こちらのページが参考になる。(注意)Cytoscapeは3.1.1以降が対応しており、Java8は未検証のようなので7が望ましい。

2.Rで実行する為に開発された関数を準備する

　RのグラフオブジェクトからCytoscape用のJSONに変換する関数などが準備されている。こちらのGitHubにあるtoCytoscape.Rがその関数。任意のフォルダに保存しておく。
　

3.シンプルなネットワーク図を描いてみる

　スタイルなどの設定を行わず簡単な図をRから描いてみる。

library(RJSONIO)
library(igraph)
library(httr)

# 基本設定
source('toCytoscape.R') # 必用な関数を読み込む
port.number = 1234 # デフォルトのポート
base.url = paste("http://localhost:", toString(port.number), "/v1", sep="") # v1はAPIのバージョン
print(base.url)

# グラフデータの準備
graph1 <- barabasi.game(200) # Barabasi-Albert model
cygraph <- toCytoscape(graph1) # Cytoscape用のJSONに変換

# Cytoscapeに送信
send2cy(cygraph=cygraph, style.name="", layout.name="force-directed")

表示された。
f:id:deta:20141227160956j:plain
動作確認というレベルだが、このように簡単に連携させることが出来る。

4.スタイルなどの設定

　任意のスタイルを適用する場合、スタイルを準備→JSON変換→APIでPOST、スタイル名を指定して描画、スタイルの適用といった流れになる感じ。次のように記述すると、CytoscapeのStyleに"MyStyle1"が登録され、そのスタイルを適用した図を描くことが出来る。

style.name = "MyStyle1"

# Defaults
def.background.color <- list(
  visualProperty = "NETWORK_BACKGROUND_PAINT",
  value = "#000000"
)

def.node.color <- list(
  visualProperty = "NODE_FILL_COLOR",
  value = "#ffffff"
)

def.edge.color <- list(
  visualProperty = "EDGE_STROKE_UNSELECTED_PAINT",
  value = "#ffffff"
)

defaults <- list(def.background.color, def.node.color, def.edge.color)

# Mappings
node.label = list(
  mappingType="passthrough",
  mappingColumn="name",
  mappingColumnType="String",
  visualProperty="NODE_LABEL"
)

mappings = list(node.label)

style <- list(title = style.name, defaults = defaults, mappings = mappings)
style.JSON <- toJSON(style)
style.url = paste(base.url, "styles", sep="/")
POST(url=style.url, body=style.JSON, encode = "json")

send2cy(cygraph=cygraph, style.name=style.name, layout.name="force-directed")

f:id:deta:20141227165643j:plain
MyStyle1が追加され、そのスタイルにあった図が表示された。

5.スタイル設定について

　スタイル設定についての詳細なドキュメントはまだない。APIを使い、スタイルの設定をGETして参考にするのが早いらしい。APIのドキュメントはこちらのページにある。APIは、おおまかにapply、networks、session、styleといった種類がある。スタイルのデフォルトは、

http://localhost:1234/v1/styles/{スタイル名}/defaults

で確認することが出来る。APIの実行や結果確認にはAdvanced REST client for Chromeが便利らしいので、これを使っている。デフォルトで登録されているスタイル名は、

http://localhost:1234/v1/styles

で確認出来る。


参考

Ruby Geocoderを使う

Ruby

　Tableauをちょいちょい使うことがあり地図表示用に緯度経度データを準備することにした。（Tableauは県名などから緯度経度情報に変換してくれる機能があるが、保有しているデータの市町村名だと結構厳しい場合がある ○船橋 ×船橋市）RubyでGeocoderってのが便利そうなので早速使ってみた。※数年前は東大CSVアドレスマッチングサービスを使用していた。

準備したRubyコード（郵便番号から住所と緯度経度を取得する）

　開発した環境はWin7、Rubyは1.9.3p448

# -*- encoding: Shift_JIS -*-
require 'geocoder'
Geocoder.configure(:language  => :ja,  :units => :km, :lookup => :google) # 設定
open("./zip.txt" ) {|fin|
  fin.each {|zip|
    res_ad = Geocoder.search(zip.chomp.encode('utf-8')) # 一応utf8にしとく
    res = [zip.chomp]
    res << res_ad[0].address_components.map{|ad| ad['long_name']}[1..-1].reverse.join(' ').encode('Shift_JIS')
    res << res_ad[0].geometry['location'].values.join("\t").encode('Shift_JIS')
    puts res.join("\t")
  }
}

　geocoderは地名やIPアドレス、緯度経度などを与えて検索することも出来る。また様々なAPIに対応している。デフォルトはGoogle Geocoding APIでGeocoder.configureで変更することができる。Google Geocoding APIの場合は1日あたりの位置情報リクエストが 2,500回に制限されているらしいので、その場合は他API用にコードを書き換える。

proxy環境の場合は、設定にproxyサーバー情報を追加する。こんな感じ。

Geocoder.configure(:language  => :ja,  :units => :km, :lookup => :google, :http_proxy => 'ユーザーID:パスワード@xxx.yyy.zzz:8080', :timeout => 10 )

取得する郵便番号のリスト（zip.txt）

103-0001
103-0002
103-0003
103-0004
103-0005

結果

103-0001	日本 東京都 中央区 日本橋小伝馬町	35.691126,139.7788326
103-0002	日本 東京都 中央区 日本橋馬喰町	35.6936397,139.7819296
103-0003	日本 東京都 中央区 日本橋横山町	35.6931621,139.7835989
103-0004	日本 東京都 中央区 東日本橋	35.6918428,139.7854238
103-0005	日本 東京都 中央区 日本橋久松町	35.6897153,139.7841379


参考

shinyIncubatorパッケージで計算中にメッセージ表示

R

　shinyIncubatorパッケージはShinyで処理中にメッセージを表示することができる。Shinyで処理に時間がかかる場合って不安になるからね。
f:id:deta:20140520054401p:plain
これはShiny - GalleryのShiny - Progress exampleで紹介されている。shinyIncubatorパッケージは正式にサポートされてないとのこと。rstudio/shiny-incubator · GitHub


shinyIncubatorパッケージはCRANにない。GitHubからインストールする。

devtools::install_github("shiny-incubator", "rstudio")

プロキシ環境の場合には、前回投稿したこちらの記事のinstall_githubの場合を参考にする。


準備したコードはShiny - Progress exampleのメッセージを日本語に書き換えたもの。

ui.R

library(shinyIncubator)

shinyUI(fluidPage(
  progressInit(),
  h1("Progress demo"),
  sidebarLayout(
    sidebarPanel(
      actionButton("go", "Run")
    ),
    mainPanel(
      plotOutput("plot")
    )
  )
))

server.R

library(shinyIncubator)

shinyServer(function(input, output, session) {
  output$plot <- renderPlot({
    if (input$go == 0)
      return()

    # Wrap the entire expensive operation with withProgress
    withProgress(session, {
      setProgress(message = "計算中。。。ちょっと待ってね", 
        detail = "計算処理に少し時間がかかります。")
      Sys.sleep(2) # 2秒待つ
      setProgress(detail = "メッセージ1")
      Sys.sleep(2) # 2秒待つ
      anotherExpensiveOperation() 
      Sys.sleep(2) # 2秒待つ
      setProgress(detail = "メッセージ2")
      Sys.sleep(2) # 2秒待つ
      plot(rnorm(100), rnorm(100))
    })
  })

  anotherExpensiveOperation <- function() {
    withProgress(session, min = 0, max = 10, {
      setProgress(message = "計算中の詳細")
      for (i in 1:10) {
        setProgress(value = i)
        if (i == 7)
          setProgress(detail = "すみません。ちょっと時間がかかっています。")
        Sys.sleep(0.3)
      }
    })
  }
})

Shinyでもd3.jsを容易に使えるd3Networkパッケージを試した

R

　最近、Rのd3Networkパッケージを使っている記事をいくつか目にした。d3Networkパッケージは、Rからd3.jsのnetwork、tree、dendogram、Sankeyを容易に描画することが出来るパッケージのようだ。great!!!　そして開発者のページを見てみると「d3Network in Shiny web apps」なるものが！　ということでさっそく試した。試したデータは普段分析している比較的小規模なネットワークのデータ。※ShinyはRで簡単にWebアプリ化することができるパッケージ。

> nd <- read.table("testdata.tsv", header=TRUE)
> head(nd)
  Source Target
1      1      2
2      1     11
3      1     21
4      2     21
5      3      4
6      3     11
>


　準備したui.Rとserver.Rは次の通り。開発者のサンプルでは不透明度合のopacityが引数になっていた。今回は斥力のchargeとノードの文字サイズのfont sizeをパタメータに追加してみた。また、shinyパッケージのfluidPage関数、fixedRow関数を使いUIを変えてみた。

ui.R

library(shiny)

shinyUI(fluidPage(

  # Load d3.js
  tags$head(tags$script(src = 'http://d3js.org/d3.v3.min.js')),

  # Application title
  titlePanel('d3Network Shiny Example'),


  # Sidebar with a slider input for node opacity
  fixedRow(
    fixedRow(
      column(3, offset=1,
        sliderInput(inputId='sld_opc', label = 'node opacity',
          min = 0, max = 1, step = 0.01, value = 0.5
        )
      ),
      column(3,
        sliderInput(inputId='sld_charge', label = 'link charge',
          min = -500, max = -10, step = 10, value = -200
        )
      ),
      column(3,
        selectInput(inputId='sel_fontsize', label = 'font size',
          choices = c(8,10,12,14,16,18,20,24),
          selected= 12
        )
      )
    ),
    hr(),
    # Show network graph
    fixedRow(
      column(12,
        htmlOutput('networkPlot')
      )
    )
  )
))


server.R

library(d3Network)

nd <- read.table("testdata.tsv", header=TRUE)

shinyServer(function(input, output) {
    output$networkPlot <- renderPrint({
      d3SimpleNetwork(nd, width=800, height=500, 
          standAlone=FALSE, opacity = input$sld_opc,
          charge = input$sld_charge, fontsize = as.numeric(input$sel_fontsize),
          parentElement = "#networkPlot")
    })
})


　実行するとブラウザではこんな感じに表示される。

library(shiny)
runApp(appDir="d:/r/ShinyApp/d3Network") # shiny実行

f:id:deta:20140513230959p:plain
　当然ShinyServerでも動作する。

Rコンソールから試したい場合

　htmlファイルに出力し、ブラウザを起動するといった流れの場合は下記のようにすれば良い。

library(d3Network)

Source <- c("A", "A", "A", "A", "B", "B", "C", "C", "D")
Target <- c("B", "C", "D", "J", "E", "F", "G", "H", "I")
NetworkData <- data.frame(Source, Target)
d3SimpleNetwork(NetworkData, width=1000, height=1000, standAlone=TRUE, file="aaa.html")
browseURL("aaa.html")

　出力したhtmlファイルを共有するってのもあり。


参考

Levenberg-Marquardtアルゴリズムで非線形回帰分析

R

　Rでminpack.lmパッケージのnls.lm関数を使うとLevenberg-Marquardt法で非線形回帰分析を行うことができる。言いかたはマルカート法、マーカート法？、フランスの方のようなのでマーカール法とも言われているようだ。私はマルカートで覚えていた。Levenberg-Marquardt法は非線形最小二乗問題を解く手法として広く使われている。最急降下法・ニュートン法を組み合わせた方法で現在の解が正解から遠い場合は遅いが収束することが保証されている最急降下法と同じように動作し、正解から近い場合はニュートン法を実行するとのこと。

非線形回帰を行うことになった背景

　最近ECサイトの分析を行うことが多く、よく言われているベキ乗則を随所に見かける。データ例

> head(d)
   X       Y
1  2 9794688
2  3 5973376
3  9 2946944
4 10 2649600
5 12 2520064
6 15 2168766
> plot(d$X, d$Y)

f:id:deta:20140508043056p:plain
これは、とあるカテゴリーのランキングと指標○○の関係。こんな関数で非線形回帰分析を行う。
y=ax^b
yは指標○○、xはランキング、aは切片(1位の指標○○値)、bは減衰項

nls.lm関数を実行する

> library(minpack.lm)
> res_lm <- lm(formula=log(Y)~log(X), data=d) # 初期パラメータ設定の為に対数変換して線形回帰
> lm_pa <- exp(res_lm$coefficients[1]) # aの初期パラメータ
> lm_pb <- res_lm$coefficients[2] # bの初期パラメータ
> x <- d$X # 説明変数
> obs <- d$Y # 目的変数
> pred <- function(parS, xx) parS$a * xx ^ parS$b # モデル式
> resid <- function(p, observed, xx) observed - pred(p,xx) # 残差
> parStart <- list(a=lm_pa, b=lm_pb) # 初期パラメータ
> nls.out <- nls.lm(par=parStart, fn=resid, observed=obs, xx=x, control=nls.lm.control(maxiter=1024,nprint=1))
It.    0, RSS = 1.59819e+15, Par. =  1.29507e+08   -1.51893
It.    1, RSS = 1.2133e+14, Par. =  3.83967e+06   -1.26897
It.    2, RSS = 9.01298e+12, Par. =  2.451e+07   -1.18445
It.    3, RSS = 4.13613e+12, Par. =  1.32097e+07  -0.719333
It.    4, RSS = 1.38756e+12, Par. =  1.64626e+07  -0.823374
It.    5, RSS = 1.34457e+12, Par. =  1.66165e+07  -0.817677
It.    6, RSS = 1.34456e+12, Par. =  1.66199e+07  -0.817936
It.    7, RSS = 1.34456e+12, Par. =  1.66198e+07   -0.81793
> nls.out
Nonlinear regression via the Levenberg-Marquardt algorithm
parameter estimates: 16619797.5617096, -0.817929785073928 
residual sum-of-squares: 1.345e+12
reason terminated: Relative error in the sum of squares is at most `ftol'.
> 
> nls_a <- nls.out$par$a
> nls_b <- nls.out$par$b
> cor(nls_a*x^nls_b, obs)
[1] 0.995052
> 
> plot(d$X, d$Y)
> x <- seq(min(d$X), max(d$X), by=1)
> x.pred <- nls_a*x^nls_b
> lines(x, x.pred, col=2)
>

f:id:deta:20140508044224p:plain
あてはまり度合も良好で説明力も高そう。今回の非線形回帰の初期値はモデル式を対数変換して線形回帰した回帰係数を使用した。Levenberg-Marquardt法で得られた回帰係数はRで一般的に知られているnls関数で得られる値とほとんど変わらない。

nlsだとエラーが出る場合が。。。

　Levenberg-Marquardt法を行う前はnls関数を使用していた。nls関数の場合は分析するデータにより下記のようにエラーが出る場合がある。algorithm="port"で改善する場合はあったが全てのエラーは解消しなかった。今のところnls.lm関数ではエラーは出てない。

Error in nls(formula = Y ~ a * X^b, data = d, start = list(a = exp(res_lm$coefficients[1]),  : 
   勾配が特異です 

Error in numericDeriv(form[[3L]], names(ind), env) : 
   モデルを評価する際，欠損値または無限大が生成されました

参考