pythonから統計解析ソフトRを利用する方法・第３回

2018/03/07

　統計解析ソフトR(以下、統計R)をpythonから利用するため
pyrcmd.py というライブラリを作りました。

　その pyrcmd の使い方について記します。

　今回は、統計処理した結果をhtmlやdocx(ワード文書)にする方法を取り上げます。

　なお、当サイトで紹介したスクリプトは
pyrcmd03.zip に入れてあります。

　スクリプトを動かすのに必要な環境については
pythonから統計解析ソフトRを利用する方法・第１回を参照してください。

《このページの目次》

1.　統計レポート生成の方針
2.　markdownのテンプレートを使った文書の作成
3.　基本的な統計レポートの作成
4.　条件分岐の例
- (1)　処理の概要
- (2)　サンプルのスクリプト

1.　統計レポート生成の方針

　統計処理した結果をhtmlとして書き出すのに、
pythonの markdown というライブラリをつかいます。

　これは、markdown形式で書かれた原稿をhtml形式に変換するためのライブラリです。

　また、docx(ワード文書)を生成するのには、
pypandoc というpythonライブラリを用います。

　これは、pandocという外部コマンドをpythonから操作するものです。

　pandocには、markdownの原稿などを様々なファイルに変換する機能があります。

　加えて、jinja2というテンプレートエンジンの利用も試します。

　統計Rの方には、統計レポートを生成するための rmarkdown という
とても強力なパッケージがあります。

　html, docx のほか、LaTeX経由で pdf を生成できます。

　markdownの原稿中に統計Rのスクリプトを書き込んでおいて、
スクリプトの実行とレポートの生成を一体的に行うものです。

　スクリプトと原稿を別々にではなく統合して管理できるのが魅力です。

　私も実践で便利に使わせてもらっています。

　ただ、実感では「条件分岐を扱いにくい」です。

　有意性が検証されたときは「○○○」、されなかった場合は「△△△」
そんな切り替えが、できないわけではありませんが分かりにくくなります。

　その点、pythonライブラリの jinja2 では条件分岐が自然に行えます。

　ということで、統計Rとpythonの両スクリプトを使いながら
レポート生成を試みます。

目次に戻る

2.　markdownのテンプレートを使った文書の作成

　markdown記法については触れないので、別のサイトを参考にして下さい。

　レポート作成の手順は、おおよそ次のようになります。

markdown原稿のテンプレート(ひな形)を用意する。
統計的な数値等は、後から入れられるように適当な名前で書いて置く。
テンプレートに入れ込む数値等をpythonスクリプトで算出または用意する。
テンプレートに数値等を入れ込み、markdown原稿を完成させる。
pythonに標準で用意されている format() あるいは jinja2 を利用。
完成したmarkdown原稿を html あるいは docx に変換する。

　上の手順の具体例を記していきます。

(1)　format() を使った数値の入れ込み

　まずは、統計処理を伴わない簡単な例です。

　下のような３行の箇条書きのテンプレートがあるとします。
(markdown原稿のテンプレートです。)

- リンゴの値段: {apple}
- バナナの値段: {banana}
- サクランボの値段: {cherry}

　apple, banana, cherry が中括弧の {……} で囲まれていますが、
そこに数値を入れれば markdownの原稿が完成です。

　原稿を完成させるための pythonスクリプトは次のとおり。

-------- ここから
# encoding: cp932
import pyrcmd as prc
md_tpl = u'''\
- リンゴの値段: {apple}
- バナナの値段: {banana}
- サクランボの値段: {cherry}
'''
md_src = md_tpl.format(apple=500, banana=450, cherry=600)
print(md_src)
-------- ここまで

　上記スクリプトの表示は下のようになります。

- リンゴの値段: 500
- バナナの値段: 450
- サクランボの値段: 600

　format(……) の引数を記述するところは、
予め辞書(hash)の形でデータを用意して置いて次のようにも書けます。
(テンプレートは同じなので省略します。)

　　（前略）
dct = {'apple': 500, 'banana': 450, 'cherry': 600}
md_src = md_tpl.format(**dct)
print(md_src)

　辞書を示す変数 dct の前に ** を置きます。

　この辞書を利用する方法は、統計処理で得た数値等を入れ込む際に
なかなか有効に機能します。

目次に戻る

(2)　jinja2を使った数値の入れ込み

　jinja2 の場合は、文字列をダイレクトにテンプレートとして用いるのではなく、
Template() を呼び出して独自の処理を施した上でテンプレートを生成します。

　そうして生成したテンプレートに対し render() を適用すれば
markdown原稿の完成です。

　format() を使う場合は {apple} のように中括弧で囲みましたが、
jinja2 の方では {{apple}} のように中括弧を二重にします。

　前述の apple, banana, cherry の例を jinja2 で処理すると次のとおり。

-------- ここから
# encoding: cp932
from jinja2 import Template
import pyrcmd as prc
md_tpl = Template(u'''\
- リンゴの値段: {{apple}}
- バナナの値段: {{banana}}
- サクランボの値段: {{cherry}}
''')
dct = {'apple': 500, 'banana': 450, 'cherry': 600}
md_src = md_tpl.render(dct)
print(md_src)
-------- ここまで

　formatのときは format(**dct) のように ** を置きましたが、
renderの場合は render(dct) とすればＯＫです。

目次に戻る

(3)　htmlの生成

　format() を使ったにしろ jinja2 を使ったにしろ、
markdownの原稿(完成版)が変数 md_src に代入されていれば
それを html に変換するのは難しくありません。

import markdown as md

　上記を pythonスクリプトの最初の方に書いておいたうえで、
下のように記述すれば htmlのbody部分を得られます。

body = md.markdown(md_src)

　これは、あくまで body部分です。

　<head>部分を含むスタンドアロンのhtmlにするには自前で加工します。

　pyrcmd.py では html_template という文字列を定義しており、
簡易な html でよければ利用できます。

　html_template には三つの %s が含まれていて、
それぞれ charset, title, body に置き換えられます。

　説明書きを並べているだけだと分かりにくいのでサンプルを掲げます。

　test.htm を utf-8 で書き出します。

-------- ここから
# encoding: cp932
import markdown as md
import pyrcmd as prc

md_tpl = u'''\
- リンゴの値段: {apple}
- バナナの値段: {banana}
- サクランボの値段: {cherry}
'''
dct = {'apple': 500, 'banana': 450, 'cherry': 600}
md_src = md_tpl.format(**dct)

title = u'果物の値段'
body = md.markdown(md_src)
html = prc.html_template % ('utf-8', title, body)
prc.write_text("test.htm", html, 'utf-8')
-------- ここまで

目次に戻る

(4)　docxファイルの生成

　markdownの原稿(テンプレートではなく完成版)が
input.md というファイルに書かれているものとします。

　docxファイル(ワード文書)を生成するためには、
import pypandoc を記述しておいたうえで下のようにします。

pypandoc.convert_file('input.md', 'docx',
    outputfile='test.docx', format='md')

　こうすると test.docx が生成されます。

　format='md' は、素材である input.md が
markdown形式であることを明示するものです。

　気をつけなければいけないのは input.md の中身が
utf-8 で書かれている必要があることです。

　pandoc というコマンドは、utf-8 以外の文字エンコーディングを扱えません。

　ちなみに、pandoc は、htmlファイルを docxファイルに変換する機能も持っています。

　下のようにします。

pypandoc.convert_file('input.htm', 'docx',
    outputfile='test.docx', format='html')

　果物の値段を示す docxファイル生成の pythonスクリプトは下のとおり。

　jinja2 を使っています。

-------- ここから
# encoding: cp932
import os
from jinja2 import Template
import pypandoc
import pyrcmd as prc

md_tpl = Template(u'''\
- リンゴの値段: {{apple}}
- バナナの値段: {{banana}}
- サクランボの値段: {{cherry}}
''')
dct = {'apple': 500, 'banana': 450, 'cherry': 600}
md_src = md_tpl.render(dct)

md_file = 'temp.md'
prc.write_text(md_file, md_src, 'utf-8')
pypandoc.convert_file(md_file, 'docx',
    outputfile='test.docx', format='md')
os.unlink(md_file)
-------- ここまで

目次に戻る

3.　基本的な統計レポートの作成

　レポート作成の基本的な手順は、これまで述べてきた事柄に即して行います。

　その際、ちょっと便利になるかも？という関数を pyrcmd.py に組み込みました。

　それらの使い方などを紹介します。

(1)　辞書データを生成するための format_od() の役割

　「第２回」で t検定の例を取り上げました。

　身長の平均値が男女間で有意に異なるか否かを見るものです。

　統計Rの側で t検定を行い、
その結果を python側で受け取るところまでを再掲します。

　python側では、rpy2のオブジェクトを
pandasのオブジェクト(Seriesなど)に変換しています。

-------- ここから
# (encoding: cp932)
import pandas as pd
from rpy2.robjects import r
from textwrap import dedent
import pyrcmd as prc
prc.set_charset('cp932')

src_file = 'data.csv'
rscript = dedent(u'''\
    dtf <- read.csv("%s", header=T, fileEncoding="%s")
    male = subset(dtf, 性別 == "男性")
    female = subset(dtf, 性別 == "女性")
    flag <- FALSE  # 不等分散を仮定
    t_test <- t.test(male$身長, female$身長, var.equal=flag)
    t_test
    ''') % (src_file, prc.get_charset(src_file))
sink_data = prc.rexec(rscript)  # Rコマンドの実行
t_test = r['t_test']  # rp2オブジェクト(ListVector)として取得
pod = prc.rp2pd(t_test)  # ListVectorを変換(OrderedDictになる)
-------- ここまで

　変数 pod は、数値、文字列または Seriesからなる OrderedDict です。

　これを材料にして markdownテンプレートに数値等を入れ込むのは
少々やっかいです。

　Series が使いにくい。

　たとえば、t検定の統計量 t を参照するのに
pod['statistic']['t'] と書かなければなりません。

　podに対して添え字を二重に指定する必要があるわけです。

　また、有意確率のp値は、pod['p.value'] として参照します。

　p.value のピリオド記号は
統計Rの記法をそのまま継承したものですが、
pythonではピリオド記号が特別の意味を持つので
変数名の一部としては使えません。

　つまり、テンプレートの中に {p.value とか {{p.value}} と書くと
エラーになるか、もしくは、意図とは別の意味になってしまいます。

　というように、変数 pod をテンプレートに適用するのは厄介です。

　そこで、podを適用しやすい形に加工するための format_od() を設けました。

fod = prc.format_od(pod)

　上のようにすると、テンプレートに対して次のように記述することで
markdownの原稿を生成できます。

md_src = tpl.format(**fod)
md_src = tpl.render(fod)  # jinja2の場合

　もちろん、テンプレートは予め作成しておく必要があります。

　fod も OrderedDict ですが、中身は下のとおり(t検定の結果)。

t: 10.3735268626
df: 386.779600635
p_value: 2.08272110497e-22
conf_int: [5.36219942931, 7.87015536867]
mean_of_x: 166.715656566
mean_of_y: 160.099479167
difference_in_means: 0.0
alternative: two.sided
method: Welch Two Sample t-test
data_name: male$身長 and female$身長

　fodを使えば、テンプレートを書く際に
{df} とか {p_value} などと記述することができます。

目次に戻る

(2)　format_od() の変換ルール

　fod = format_od(pod) が pod をどのように変換するのか、
そのルールを列記します。

OrderedDictであるpodのkeyについて、ピリオドやスペースなど
pythonの変数名として使えない文字をアンダーラインに置き換える。
p.value → p_value　｜　mean of x → mean_of_x
添え字が二重になるものは、一重になるよう変換。
pod['a']['t'] → fod['t'] ('a' は無視してしまう)
後の方に pod['b']['t'] があれば、それは fod['b_t'] に変換。
prefixを指定したときは、keyの先頭に指定文字をつける。
fod = format_od(pod, prefix='T') とすれば、
fod['Tdf'] とか fod['Tp_value'] などとなる。
DataFrameは markdownのパイプテーブルに変換する。
podがOrderedDictの入れ子構造になっている場合は、子dictも変換する。
なお、先の「添え字が二重のものを一重にする」は、親と子各々に適用し、ごっちゃにはしない。

　prefixを指定できるようにしてあるのは、
複数の検定結果を扱えるようにするためです。

　F検定の結果が変数 f_fod に入っていて、
t検定の結果が t_fod に入っているとします。

　この二つの OrderedDict を一つのテンプレートに適用しようとすると、
まず dict の連結を行うことになります。

　しかし、どちらにも p_value があるので、連結すると一方のデータが失われます。

　そこで、接頭辞 prefix を指定して

f_fod = format_od(f_pod, prefix='F')
t_fod = format_od(t_pod, prefix='T')

　上のようにすれば
f_fod['Fp_value'] および t_fod['Tp_value'] となって
連結しても差し支えないことになります。

　すべての key に接頭辞がつくのは煩わしい感じもしますが、
テンプレートを一つにまとめて管理したい場合は、このように prefix を使えます。

　ちなみに、辞書の連結は update() で行うことができます。

d1 = {'a': 1}
d2 = {'b': 2}
d1.update(d2)
print(d1)  # {'a': 1, 'b': 2}

　format_od() が DataFrame をどのように扱うかは、後の方で述べます。

目次に戻る

(3)　t検定の結果をhtmlファイルとして書き出す

　ここでは markdownのテンプレートを format() でしょりします(jinja2ではない)。

　統計Rの出力結果 sink_data を html の最後のところに挿入します。

　ちょっと長いですが、省略せずに pythonスクリプト全体を掲げます。

 1# prc01s.py (encoding: cp932)
 2import pandas as pd
 3from rpy2.robjects import r
 4from textwrap import dedent
 5import markdown as md
 6import pyrcmd as prc
 7prc.set_charset('cp932')
 8
 9src_file = 'data.csv'
10rscript = dedent(u'''\
11    dtf <- read.csv("%s", header=T, fileEncoding="%s")
12    male = subset(dtf, 性別 == "男性")
13    female = subset(dtf, 性別 == "女性")
14    flag <- FALSE  # 不等分散を仮定
15    t_test <- t.test(male$身長, female$身長, var.equal=flag)
16    t_test
17    ''') % (src_file, prc.get_charset(src_file))
18md_tpl = dedent(u'''\
19    # 身長の平均値に関する検証\n
20    　身長の平均値が男女間で異なるかを見るため t検定(両側検定)を行った。\n
21    　主な数値は次のとおり。\n
22    - 統計量 t値: {t}
23    - 自由度(df): {df}
24    - 有意確率(p値): {p_value}
25    - 平均値: 男性 {mean_of_x:.1f}, 女性 {mean_of_y:.1f}
26    - 95％信頼区間: {conf_int[0]:.4f} -- {conf_int[1]:.4f}\n
27    　有意確率が低いことから(p<0.05)、平均値が等しくないと判断される。\n
28    　なお、統計解析ソフトRの出力結果は下のとおり。\n
29    ~~~~
30    {sink}
31    ~~~~\n
32    ----\n\n
33    ''')
34
35sink_data = prc.rexec(rscript)  # Rコマンドの実行
36t_test = r['t_test']  # rp2オブジェクト(ListVector)として取得
37pod = prc.rp2pd(t_test)  # ListVectorを変換(OrderedDictになる)
38fod = prc.format_od(pod)  # テンプレート用に変換
39fod['sink'] = sink_data  # 統計Rの出力を組み入れる
40md_src = md_tpl.format(**fod)  # markdown原稿の生成
41title = u't検定'  # htmlのタイトル
42body = md.markdown(md_src, ["fenced_code"])
43html = prc.html_template % ('utf-8', title, body)
44prc.write_text("prc01s.htm", html, 'utf-8')

　markdownのテンプレートの中で男女それぞれの平均値を示す箇所があります。

男性 {mean_of_x:.1f}, 女性 {mean_of_y:.1f}

　この :.1f は表示の書式を指定するもので、小数点以下の桁数指定です。

　単に {mean_of_x} とすると、ありったけの桁数が表示されてしまうので
小数点以下１桁までの表示にするためこのように記述しています。

　統計Rの出力結果 sink_data は、プログラムソースコードと同じように
<pre><code> ………… </code></pre> で囲む形にするため
markdown原稿の中では ~~~~~~~~ という行で囲んでいます。

　このような囲まれたコードを処理する場合は、

body = md.markdown(md_src, ["fenced_code"])

　上のように ["fenced_code"] というオプションを指定します。

　なお、sink_data は、podを変換した元々の fod には含まれていないので
fod['sink'] = sink_data のようにして組み入れます。

　そうすることで、テンプレートに {sink} と書けるようになります。

目次に戻る

(4)　DataFrameのパイプテーブルへの変換

　「第２回」でカイ２乗検定を取り上げました。

　「身長区分」×「性別」のクロス集計表を検定するものです。

　検定結果にはカイ２乗統計量、自由度、p値などのほかに、
観測度数や調整済み残差などの表データ(テーブル)も含まれています。

　表データをpython側で受け取ると DataFrame になりますが
これを markdownの原稿に組み込む一つの方法は
パイプテーブルに変換する、というものです。

　セルとセルを半角の縦線記号(パイプ記号)で区切る形式です。

|品名|価格|
|-----|----:|
|リンゴ|500|
|バナナ|450|
|サクランボ|600|

　必ずしも、パイプ記号の位置をそろえる必要はありません。

　２行目の |-----|----:| はヘッダと本体の区切りで、
コロン記号の位置により左揃え・中央揃え・右揃えを指定できます。

　上記は、１列目は位置の指定なしですが、２列目は右揃えを指定しています。

　変数 pod(OrderedDict型)の要素に DataFrame があると、
fod = format_od(pod) と変換した場合、
DataFrame はパイプテーブルに変換されます。

　DataFrameの中の数値データは右揃え、それ以外は指定なしです。

　「身長区分」×「性別」の場合でいうと、
観測度数の fod['observed'] が下のようになります。

||男性|女性|
|-----|----:|----:|
|155未満|6|40|
|155～165|72|109|
|165～175|100|41|
|175以上|20|2|

　パイプテーブルへの変換についてはオプション指定ができますが、
くどい説明は後回しにして、カイ２乗検定の結果をレポート化してみます。

　jinja2 を使って markdownの原稿を生成し、
pypandocによって docxファイルを書き出します。

　長いですが全部を掲げます。

 1# prc02s.py (encoding: cp932)
 2import os
 3import pandas as pd
 4from rpy2.robjects import r
 5from jinja2 import Template
 6from textwrap import dedent
 7import pypandoc
 8import pyrcmd as prc
 9prc.set_charset('cp932')
10
11src_file = 'data.csv'
12rscript = dedent(u'''\
13    dtf <- read.csv("%s", header=T, fileEncoding="%s")
14    dtf$性別 <- factor(dtf$性別, levels=c("男性", "女性"),
15      labels=c("男性", "女性"))
16    cc <- cut(dtf$身長, breaks=c(-Inf, 155, 165, 175, Inf),
17        labels=c("155未満", "155～165", "165～175", "175以上"), right=FALSE)
18    dtf <- cbind(dtf, 身長区分=cc)
19    tbl <- table(dtf$身長区分, dtf$性別)
20    chi_test <- chisq.test(tbl)
21    chi_test
22    ''') % (src_file, prc.get_charset(src_file))
23md_tpl = Template(dedent(u'''\
24    # 「身長区分」×「性別」のカイ２乗検定結果\n
25    　男女間の身長の違いをみるため、下の観測度数の表をカイ２乗検定にかけた。\n
26    {{observed}}\n
27    　検定結果の主な数値は次のとおり。\n
28    - 統計量(カイ２乗値): {{X_squared}}
29    - 自由度(df): {{df}}
30    - 有意確率(p値): {{p_value}}\n
31    　有意確率の値が低いことから(p<0.05)、
32    「身長区分」と「性別」の関係は独立なものではなく、
33    有意な関連があると判断される。\n
34    　そこで、調整済み残差をみると下のとおり。\n
35    {{stdres}}\n
36    　これをみると、165未満では女性が有意に多く、
37    一方、165以上では男性が有意に多いことが分かる。\n
38    '''))
39
40prc.rexec(rscript, sink=False)  # Rコマンドの実行
41chi_test = r['chi_test']  # rp2オブジェクト(ListVector)として取得
42pod = prc.rp2pd(chi_test)  # ListVectorを変換(OrderedDictになる)
43fod = prc.format_od(pod)
44md_src = md_tpl.render(fod)  # markdown原稿の生成
45md_file = 'temp.md'
46prc.write_text(md_file, md_src, 'utf-8')  # markdown原稿をファイルに
47pypandoc.convert_file(md_file, 'docx',
48    outputfile='prc02s.docx', format='md')
49os.unlink(md_file)

　pythonのmarkdownライブラリを用いて htmlにする場合は
import markdown as md を前の方に書いておいた上で下のようにします。

title = u'カイ２乗検定'  # htmlのタイトル
body = md.markdown(md_src, extensions=["tables"])
html = prc.html_template % ('utf-8', title, body)
prc.write_text("prc02s.htm", html, 'utf-8')

　md.markdown(……) では tables というオプションを指定しないと、
パイプテーブルがうまく処理されないようです。

目次に戻る

(5)　パイプテーブルに変換する際のオプション

　pyrcmd.py の中で md_table() という関数を定義しています。

　一つの pandas.DataFrame をパイプテーブル(文字列)に変換するものです。

md_table(dtf, align=None, index=True)

　オプションに align, index があります。

　index=False とすれば、DataFrameの行ラベルを取り除いた上で
パイプテーブルに変換します。

　align は、各セル内の表示位置を指定します。

　align = 'lcr' とすれば、第１列が左揃え、
第２列が中央揃え、第３列が右揃えになります。

　c, l, r 以外の文字(たとえば !)を書くと、
該当の列は「指定なし」となります。

　align = 'llrr!c' のように６個の文字を並べても
実際の列数が４個しかなければ、最後の方の２個は無視されます。

　align = 'lr' と２個だけ指定したときに、
実際の列数が４個なら、最後の方の２個は「指定なし」になります。

　format_od() は、DataFrameに出会ったときに
内部で md_table() を呼び出して処理します。

　format_od() でも align, index のオプションを指定できますが、
どの DataFrame にも同じ align, index が適用されるので注意してください。

　カイ２乗検定の結果が変数 pod に代入されている場合、
観測度数、期待値、残差、調整済み残差の４つの DataFrame が入っています。

fod = prc.format_od(pod, align='!ll', index=True)

　上のようにすると、４つの DataFrame に対して
同じ align, index が適用されます。

　同じでは困るという場合は、個々の DataFrame ごとに
md_table() を適用することになるでしょうか。

［補足］　format_od() の table_typeオプション

　　　　話がややこしくなって申し訳ありませんが、
実は format_od() では table_type というオプションも指定できます。

　format_od(pod, table_type='html') とすれば、
DataFrameをパイプテーブルではなく htmlに変換します。

　つまり、<table border="1"> …… </table> に変換します。

　これは pandas.DataFrame の to_html() を呼び出して変換するものです。

　format_od() において指定されたオプションのうち、
prefix, table_type 以外は to_html() にそのまま引き渡されます。

　table_type オプションのデフォルト値は 'markdown' です。

　table_type='' とすれば、format_od() は、
DataFrameに対して変換処理を行いません。
つまり、DataFrameのままになります。

［参考］　Seriesを表形式で表示したい場合

　pandasの Series をテーブルとして表示したいこともあるとおもいます。

　format_od() は、Seriesを各要素ごとに分解した上で記録するので、
一つのテーブルとして扱うことはできなくなります。

　変数 ser に Seriesオブジェクトが代入されている場合でいうと、
次のような２行でパイプテーブルに変換できます。

tbl = pd.concat([ser], axis=1).T  # DataFrameに変換
md_tbl = prc.md_table(tbl, index=False)

　統計Rの summary() をpython側で受け取ると Series になりますが、
それに上記の２行を適用すると、たとえば下のような結果が得られます。

|Min.|1st Qu.|Median|Mean|3rd Qu.|Max.|NA's|
|----:|----:|----:|----:|----:|----:|----:|
|144.8|158.67|163.6|163.46|168.02|185.2|8.0|

　要するに、SeriesをDataFrameに変換すれば md_table() で変換できます。

　Series→DataFrame の変換方法は上に掲げたものだけではないとおもいますが、
一つの例として参考にして下さい。

目次に戻る

4.　条件分岐の例

　jinja2を使った条件分岐の例を取り上げてみます。

(1)　処理の概要

　htmlファイルを生成します。

　有意性が検証されたら「○○○」の検定を行い、
されない場合は「△△△」の検定に進む、
といった実践的な例だとかなり長くなるので
「検証された、されない」のメッセージを切り替える程度にします。

　data.csv には約400人分の身長と体重のデータが入っていますが、
身長と体重それぞれについて次の処理を行います。

summaryによる分布状況の確認
ヒストグラムの作成(pngファイル)
shapiro.test による正規性の検証

　統計Rの側では、変数 xx に情報一式を記録します。

　xx はリストで、値を参照する例でいうと下のようになります。

xx$height$summary,  xx$weight$summary
xx$height$shapiro$p.value,  xx$weight$shapiro$p.value

　python側ではそれを変数 pod や fod に格納します。

fod['height']['summary']
fod['weight']['summary']
fod['height']['shapiro']['p_value']
fod['weight']['shapiro']['p_value']

　上のような入れ子構造の辞書データとして格納します。

　jinja2のテンプレートで入れ子構造の辞書データを参照するときは
{{shapiro.p_value}} のようにピリオド記号でつなぐようです。

　なお、テンプレートは一つだけ用意して、
身長と体重の各々に適用します。

　つまり、一つのテンプレートを二度使います。

目次に戻る

(2)　サンプルのスクリプト

　統計Rのスクリプトと markdownのテンプレートは、
ほんとは別ファイルにしておくのがいいのだとおもいますが、
ばらばらだとサンプルとして把握しにくいので一体化してあります。

　70行くらいの長いスクリプトですが、全部を掲げます。

 1# prc03s.py (encoding: cp932)
 2import pandas as pd
 3from rpy2.robjects import r
 4from jinja2 import Template
 5from textwrap import dedent
 6import markdown as md
 7import pyrcmd as prc
 8prc.set_charset('cp932')
 9
10src_file = 'data.csv'
11rscript = dedent(u'''\
12    dtf <- read.csv("%s", header=T, na.string="", fileEncoding="%s")
13    jp_names <- c("身長", "体重")
14    eng_names <- c("height", "weight")
15    xx = list()  # 統計処理結果を入れるための空のリスト
16    for (i in 1:length(jp_names)) {
17        jp <- jp_names[i]
18        eng <- eng_names[i]
19        vct <- dtf[[jp]]  # 着目の列をvctに
20        vct <- vct[!is.na(vct)]  # 欠損値の除去
21        smr <- summary(vct)  # 分布の要約を得る
22        shp <- shapiro.test(x=vct)  # 正規性の検証
23        xx[[eng]] <- list(summary=smr, shapiro=shp)
24        xx[[eng]]$name <- jp
25        fname <- sprintf("prc03s_%%s.png", eng)
26        png(fname, width=240, height=240)
27        hist(vct, main=eng)
28        dev.off()  # 描画デバイスを閉じる
29        xx[[eng]]$png_name <- fname
30    }
31    ''') % (src_file, prc.get_charset(src_file))
32
33md_tpl = Template(dedent(u'''\
34    # {{name}}の分布\n
35    　{{name}}の分布状況(summary)は下のとおり。\n
36    {{summary}}\n
37    　ヒストグラムを描くと下のようになった。\n
38    ![{{name}}のヒストグラム]({{png_name}})\n
39    　正規性を検証するため shapiro.test を行ったところ、
40    次のような結果だった。\n
41    - 統計量 W: {{shapiro.W}}
42    - 有意確率 p: {{shapiro.p_value}}\n
43    {% if shapiro.p_value < 0.05 %}
44    {{name}}については、p<0.05 であることから
45    「正規分布に即している」(帰無仮説)が棄却される。\n
46    つまり、正規分布ではないと判断できる。
47    {% else %}
48    {{name}}については、p>=0.05 であることから
49    「正規分布に即している」(帰無仮説)を棄却できない。\n
50    つまり、正規分布に即しているとみることができる。
51    {% endif %}
52    --------\n\n
53    '''))
54
55prc.rexec(rscript, sink=False)  # Rコマンドの実行
56xx = r['xx']  # rpy2のオブジェクトを得る
57pod = prc.rp2pd(xx)  # pandas型などに変換
58for key in ['height', 'weight']:
59    ser = pod[key]['summary']  # summaryの結果(Series)
60    ser = ser.round(2)  # 数値のまるめ処理
61    df = pd.concat([ser], axis=1).T  # DataFrameに変換
62    pod[key]['summary'] = df
63fod = prc.format_od(pod, index=False)  # テンプレートで扱えるように変換
64
65md_src = ''
66for key in ['height', 'weight']:
67    md_src = md_src + md_tpl.render(fod[key])
68title = u'分布状況の確認と正規性の検証'  # htmlのタイトル
69body = md.markdown(md_src, extensions=["tables"])
70html = prc.html_template % ('utf-8', title, body)
71prc.write_text("prc03s.htm", html, 'utf-8')

　fod = prc.format_od(pod, index=False) の index=False は、
パイプテーブルに変換する際に行ラベルを表に入れないための指定です。

　markdownのテンプレート md_tpl は、
md_tpl.render(fod['height']) および
md_tpl.render(fod['weight']) として二度用いています。

　上に掲げたスクリプト程度であれば、jinja2 を使わず
標準の format() でも同じ処理が可能です。

　条件分岐に関係ないテンプレート部分を md_tpl1 とし、
分岐させる部分は md_tpl2 とします。

　md_tpl2 は辞書データで、md_tpl2[True], md_tpl2[False] の二つを用意します。

　あと、jinja2 において {{shapiro.p_value}} と記述したところは
{shapiro[p_value]} と書きます。

　p_value を引用符で囲んで {shapiro['p_value']} とすると
エラーになるので注意して下さい。

　zip圧縮ファイルには format() で処理するスクリプトも入っています。
(prc04s.py, prc04u.py)

　pyrcmd.py について３回にわたって書きましたが、
仕組みや機能を一通り説明できたとおもいます。

　markdown, jinja2, pandoc についてはごく単純な使い方しか紹介できず
書き足りないことが山積ですが、これで終了とします。

～　以上　～

pythonから統計解析ソフトRを利用する方法・第３回

《このページの目次》

1. 統計レポート生成の方針

2. markdownのテンプレートを使った文書の作成

(1) format() を使った数値の入れ込み

(2) jinja2を使った数値の入れ込み

(3) htmlの生成

(4) docxファイルの生成

3. 基本的な統計レポートの作成

(1) 辞書データを生成するための format_od() の役割

(2) format_od() の変換ルール

(3) t検定の結果をhtmlファイルとして書き出す

(4) DataFrameのパイプテーブルへの変換

(5) パイプテーブルに変換する際のオプション

［補足］ format_od() の table_typeオプション

［参考］ Seriesを表形式で表示したい場合

4. 条件分岐の例

(1) 処理の概要

(2) サンプルのスクリプト

1.　統計レポート生成の方針

2.　markdownのテンプレートを使った文書の作成

(1)　format() を使った数値の入れ込み

(2)　jinja2を使った数値の入れ込み

(3)　htmlの生成

(4)　docxファイルの生成

3.　基本的な統計レポートの作成

(1)　辞書データを生成するための format_od() の役割

(2)　format_od() の変換ルール

(3)　t検定の結果をhtmlファイルとして書き出す

(4)　DataFrameのパイプテーブルへの変換

(5)　パイプテーブルに変換する際のオプション

［補足］　format_od() の table_typeオプション

［参考］　Seriesを表形式で表示したい場合

4.　条件分岐の例

(1)　処理の概要

(2)　サンプルのスクリプト