RとQGISとPostgreSQLと

FOSS4G Advent Calendar 2017 18日目の記事です。

 

毎年皆さま新しい技術記事のなか、今さらながらRを始めたのでQGISプラグインから使ってみた記事です。

 

R:オープンソースの統計解析向けプログラミング言語。初心者です。

QGISオープンソースGISソフト。お世話になってます。

PostgreSQLオープンソースのデータベースソフト。お世話になってます。

Python:いろいろ使えるスクリプト言語。チューニング次第で早くなったり遅くなったり。

 

で、この方々の関係がどうなっているかとゆうと

 

R - QGIS

プロセッシングでRスクリプトを利用できる。

https://docs.qgis.org/2.18/ja/docs/training_manual/processing/r_intro.html

 

R- PostgreSQL

RPostgreSQLパッケージで利用できる。

 

R- Python

今はRを使わないでPythonで統計をやるのが流行りらしい。

Python からRを利用できるライブラリのPypeRとゆうものがあるらしい。

 

QGIS - PostgreSQL

レイヤとしてテーブル内容を読み込めるし、プラグインからSQLを実行できる。

 

QGIS - Python

PythonQGISプラグインを作成できる。

 

Python - PostgreSQL

ドライバモジュール使用で接続できる。

 

らしい。

とゆうことでQGISプラグインPythonで作成して、PyperでRを実行。R内でRPostgreSQLを使ってデータベースを更新してQGISに更新結果を反映してみる試み。

 

QGISとRとPostgreSQLはインストール済み。

 

データは国土数値情報の市区町村データと統計局estatの統計データから出生数・人口総数・婚姻数等をPostgreSQLに突っ込む。

 

f:id:norimakibros:20171208142330p:plain

 

DB 名:rsample

市区町村名テーブル:shikuchouson。これのbornに出生数、born2に推定出生数、born_diffに差異を入れる。

 

統計データテーブル:tbl2017_a。

a1101 :総人口、a4101:出生数、a9101:婚姻件数

 

これをRで出生数を目的変数、人口総数、婚姻数を説明変数とし、重回帰分析して出生数の実数、推計数、推計差分をQGISプラグインから実行してテーブルを更新して表示してみる。

 

まずQGISプラグインを作成する。作り方は矢崎さんの記事などを参考。

QGISプラグインの作り方(パッケージ生成編) - Qiita

私の環境(QGIS2.18)だとresources_rc.pyをresources.pyにリネームしないとプラグインが壊れていますとか言われた。GUIはこんな感じ。この推計出生数ボタンを押された時にRスクリプトを実行してDBを更新させる。

f:id:norimakibros:20171214154424p:plain

 

 

 次にQGIS環境変数を確認・設定。

PYTHONHOME: QGISで使用するpythonのインストール先

PYTHONPATH:pythonのsite-packagesのインストール先

R_HOME : Rのインストール先

PATH:R_HOMEを追加

 

そしてPYTHONHOMEのpythonにRのライブラリPypeRをpipでインストールする。対応するpythonのsite-packagesにインストールされるはず。

pip install pyper

インストールが正常にいき、パスがちゃんと通っていればQGISプラグイン>pythonコンソール上でimportできるはず。エラーがでたら失敗です。pythonのバージョンいろいろ入れている人はパス関係を確認してください。私はここで、はまりました。

 

f:id:norimakibros:20171214151610p:plain

 

最後にRにPostgreSQLに接続するためのRPostgreSQLパッケージをインストール。

Rを起動して

install.packages("DBI")

install.packages("RPostgreSQL")

f:id:norimakibros:20171214162755p:plain

 これで使用準備は終わり。後はRのスクリプトを作成する。今回のスクリプトは以下のようなものを born.Rとして保存。

 

#RPostgreSQLライブラリ呼び出し
require("DBI")
require("RPostgreSQL")

# DB連結
con <- dbConnect(PostgreSQL(), host="localhost", port=5434, dbname="rsample", user="postgres", password="postgres")

# データセット作成
dataset <- dbGetQuery(con,"SELECT * FROM tbl2017_a")

# lmで重回帰分析
ans <- lm(a4101~a1101+a9101, dataset)

# パラメータ取得
intercept <- as.character(ans$coefficients["(Intercept)"])
a1101 <- as.character(ans$coefficients["a1101"])
a9101 <- as.character(ans$coefficients["a9101"])

# 推計式作成
formura = paste(intercept,"+",a1101,"*a1101 +",a9101,"*a9101")

# 推計値を入力するSQL
sql = paste("UPDATE shikuchouson SET born2=a.born2 FROM",
"(SELECT scode,round(", formura, ",0) AS born2 FROM tbl2017_a) a",
"WHERE a.scode=TO_NUMBER(n03_007, '99999');")

# 実際値と推計値の差異を入力するSQL
sql2 = "UPDATE shikuchouson SET born_diff=abs(born-born2);"

# SQL実行
dbGetQuery(con,sql)
dbGetQuery(con,sql2)

 

そしてプラグイン側でこのスクリプトを読み込み実行

#---------------------------------------------------------------------------------------
# R推計値セット
#---------------------------------------------------------------------------------------

import pyper as py

def setR(self):
    success = True
    try:
        r = py.R()
        r("source(file='C:/tmp/R/QGIS/born.R')")
        QtGui.QMessageBox.information(self.dlg, u'確認', u'R推測値セットしました',         QtGui.QMessageBox.Close)
    except:
        QtGui.QMessageBox.critical(self.dlg, u'エラー', u'R推測値セットに失敗しました',         QtGui.QMessageBox.Close)
        success = False

 

import pyper as py でpyperをインポート。pyを別名にする。

r = py.R() でRのインスタンスを作成

r("source(file={ファイル名})") でRのスクリプトを実行する。

 

これはR側で全部処理した場合ですが、python側で処理したものをRとやりとりする場合は

r.assign('data', test) # rにtestをdataとして渡す

res = r.get("result") # rのresultをresとして受け取る

のような形で受け渡しできます。

 

f:id:norimakibros:20171214154714p:plain

推計値で色分け。

f:id:norimakibros:20171214154739p:plain

推計値差異で色分け。濃いとこが差異多い。

 

ちなみに北海道では札幌市を各区別としたら、旭川市が一番出生数が多いらしい。

Rは地図用のパッケージいろいろあるみたいですね。来年勉強します。