Up

チョコレート選好の一貫性

 

星野ら(2020)のTable 4によれば、「ノーズクリップなし」条件でチョコレートの選択が「ノーズクリップあり」条件でも同じである対象者は表1に示す人数であった。

表1 チョコレートの選好データ

 

ノーズクリップなし

ノーズクリップありでも

同じチョコレートを選択

ガーナ

55

31

ベネズエラ

24

8

エクアドル

16

9

 

「ノーズクリップなし」条件で選んだチョコレートが、さらに「ノーズクリップあり」条件でも選ばれる頻度は2項分布に従うと考えることができる。すなわち、「ノーズクリップなし」条件であるチョコレートが人に選ばれ、その人の内人が同じチョコレートを「ノーズクリップあり」条件で選んでいるとき、

である。ここで、はチョコレートに対して「ノーズクリップなし」で選ばれたときに「ノーズクリップあり」でも選ばれる確率である。

これらの選択確率は実験で用いられたチョコレートという共通点があり、何らかのまとまりがあると考えられる。いま、選択確率をロジット変換した値が共通の正規分布に従っていると考える(階層モデル;hierarchical model)。ロジット変換するのは、確率の変域が0から1の値であり、これらをロジット変換すると値の変域はー∞から+∞となり、正規分布が使えるからである。

ここで、

である。

式(1)の右辺の正規分布のパラメータhyperparameterと呼ばれている。これらのhyperparameterの事前分布をできるだけ平らな事前分布として次式で与える。

これらhyperparameterの事前分布は、hyperprior分布と呼ばれている。出来るだけ平らであると同時に事後分布の安定性を考えて設定した。

上の設定の確率モデルのStanスクリプトをリスト1に示す。

 

リスト1 階層モデルのStanスクリプト(persistence_hi.stan

data {

    int Ng;

    int kg;

    int Nv;

    int kv;

    int Ne;

    int ke;

}

parameters {

    real hm;

    real<lower = 0.0> hsgm;

    real hg;

    real hv;

    real he;

}

transformed parameters {

    real<lower = 0, upper = 1> pg;

    real<lower = 0, upper = 1> pv;

    real<lower = 0, upper = 1> pe;

    pg = inv_logit(hg);

    pv = inv_logit(hv);

    pe = inv_logit(he);

}

model {

    hm ~ normal(0, 10);

    hsgm ~ uniform(0, 20); 

    hg ~ normal(hm, hsgm);

    hv ~ normal(hm, hsgm);

    he ~ normal(hm, hsgm);

    kg ~ binomial(Ng, pg);

    kv ~ binomial(Nv, pv);

    ke ~ binomial(Ne, pe);

}

 

リスト1のStanスクリプトを用いて表1のデータを分析するPythonスクリプトをリスト2に示す。

Stanについては、岡本(2019で説明している。

これらのスクリプトファイルは圧縮ファイルfileschocoPersist.zipにまとめた。ダウンロード・展開すれば、自由に利用できる。

 

 

リスト2 リスト1のStanスクリプトを用いるPythonスクリプト(AnaPersist.py

import numpy as np

import matplotlib.pyplot as plt

import seaborn as sb

import pystan

import platform

 

sm = pystan.StanModel(file = 'persistence_hi.stan')

print('Compilig ended.')

 

Data = {'Ng':55, 'kg':31, 'Nv':24, 'kv':8, 'Ne':16, 'ke':9}

 

fit = sm.sampling(data = Data,

                  n_jobs = 1 if platform.system() == 'Windows' else -1)

 

print(fit)

 

pg_hat = np.median(fit['pg'])

pv_hat = np.median(fit['pv'])

pe_hat = np.median(fit['pe'])

sb.kdeplot(fit['pg'], label = 'Ghana')

sb.kdeplot(fit['pv'], label = 'Venezuela')

sb.kdeplot(fit['pe'], label = 'Equador')

plt.title('pv_hat={0:.3f},  pv_hat={1:.3f},  pe_hat={2:.3f}'.format(pg_hat, pv_hat, pe_hat))

plt.legend()

plt.savefig('Fig_persist_hi.png')

plt.show()

 

sb.kdeplot(fit['hg'], label = 'hg:Ghana')

sb.kdeplot(fit['hv'], label = 'hv:Venezuala')

sb.kdeplot(fit['he'], label = 'he:Equador')

plt.legend()

plt.savefig('FigH.png')

plt.show()

 

sb.kdeplot(fit['hm'], label = 'Hyper Mu')

plt.legend()

plt.savefig('FigHM.png')

plt.show()

 

sb.kdeplot(fit['hsgm'], label = 'Hyper Sigma')

plt.legend()

plt.savefig('FigHS.png')

plt.show()

 

 

plt.hist(fit['hsgm'], label = 'Hyper Sigma')

plt.legend()

plt.savefig('FigHSgmHist.png')

plt.show()

 

 

リスト2のスクリプトを実行すると以下の結果を得る。

図1のグラフは、3つのチョコレート、ガーナ、ベネズエラ、エクアドル、における2項分布の確率pgpvpeの事後分布である。ベネズエラの確率分布が他の2つより小さい、すなわち「ノーズクリップなし」条件から「ノーズクリップあり」条件に変わると、他のチョコレートを選びやすいことが判る。

図1

 

図2は、pgpvpeのロジット変換(式(2))の値hghvheの事後分布である。

図2

 

図3は、ハイパーパラメータの事後分布である。

図3

 

図4は、ハイパーパラメータの事後分布である。

図4

 

図4のsigmaの分布が、負の領域まで延びているのは、KDEを用いたグラフであることによる。

ヒストグラムは図5のようになる。

図5

 

分布は0以上の領域になっている。

 

 

参考文献

星野亜由美、石川華衣、岡田綾、飯田文子(2020).カカオ豆の香がチョコレート嗜好に及ぼす影響.日本官能評価学会誌、2429-42

 

 

Up