FOSS4G Advent Calendar 2017 18日目の記事です。
毎年皆さま新しい技術記事のなか、今さらながらＲを始めたのでQGISのプラグインから使ってみた記事です。
R：オープンソースの統計解析向けプログラミング言語。初心者です。
QGIS：オープンソースのGISソフト。お世話になってます。
PostgreSQL：オープンソースのデータベースソフト。お世話になってます。
Python：いろいろ使えるスクリプト言語。チューニング次第で早くなったり遅くなったり。
で、この方々の関係がどうなっているかとゆうと
R - QGIS
プロセッシングでＲスクリプトを利用できる。
https://docs.qgis.org/2.18/ja/docs/training_manual/processing/r_intro.html
R- PostgreSQL
RPostgreSQLパッケージで利用できる。
R- Python
今はRを使わないでPythonで統計をやるのが流行りらしい。
Python からRを利用できるライブラリのPypeRとゆうものがあるらしい。
レイヤとしてテーブル内容を読み込めるし、プラグインからSQLを実行できる。
ドライバモジュール使用で接続できる。
らしい。
とゆうことでQGISのプラグインをPythonで作成して、PyperでＲを実行。R内でRPostgreSQLを使ってデータベースを更新してQGISに更新結果を反映してみる試み。
QGISとＲとPostgreSQLはインストール済み。
データは国土数値情報の市区町村データと統計局estatの統計データから出生数・人口総数・婚姻数等をPostgreSQLに突っ込む。
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はこんな感じ。この推計出生数ボタンを押された時にＲスクリプトを実行してＤＢを更新させる。
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のバージョンいろいろ入れている人はパス関係を確認してください。私はここで、はまりました。
最後にRにPostgreSQLに接続するためのRPostgreSQLパッケージをインストール。
Rを起動して
install.packages("DBI")
install.packages("RPostgreSQL")
これで使用準備は終わり。後はＲのスクリプトを作成する。今回のスクリプトは以下のようなものを 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)
#---------------------------------------------------------------------------------------
# Ｒ推計値セット
#---------------------------------------------------------------------------------------
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'Ｒ推測値セットしました', QtGui.QMessageBox.Close)
except:
QtGui.QMessageBox.critical(self.dlg, u'エラー', u'Ｒ推測値セットに失敗しました', 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として受け取る
のような形で受け渡しできます。
推計値で色分け。
推計値差異で色分け。濃いとこが差異多い。
ちなみに北海道では札幌市を各区別としたら、旭川市が一番出生数が多いらしい。
Ｒは地図用のパッケージいろいろあるみたいですね。来年勉強します。
