Up

甘味度と温度の関係

糖類間の比較のダミー変数による分析

 

糖の甘味度と温度の関係が前橋（2011）に報告されている。前橋（2011）の第２図からフラクトースとグルコースの温度と甘味度の関係を読み取ったものを表１に示す。

表１　前橋（2011）の第2図からの読み取り

フラクトース

グルコース

温度℃

甘味度

温度℃

甘味度

5

1.43

3

0.59

18

1.27

19

0.57

40

1

40

0.56

60

0.79

60

0.56

 

甘味度を従属変数、温度を独立変数とする回帰モデルをフラクトースおよびグルコースに対して以下のように設定することができる。

いま、ダミー変数をフラクトースに対して１、グルコースに対して０をとるものを用意する。

である。

このとき、回帰モデルを次のようにおく。

この回帰式は、フラクトースのときは

となり、グルコースのときは

となる。

このダミー変数を用いた回帰モデルのStanスクリプトをリスト１のように用意した。

スクリプトファイルおよびデータファイルは圧縮ファイルfilesSweetTemp.zipにまとめた。ダウンロード・展開すれば自由に利用できる。

Stanの説明は、岡本（2019）で行っている。

 

リスト１　ダミー変数を用いた回帰モデルの分析Stanスクリプト（ModelDummy.stan）

data {

    int N;

    real Y[N];

    real X[N];

    real G[N];

}

parameters {

    real a1;

    real a2;

    real b1;

    real b2;

    real<lower = 0.0> sgm;

}

transformed parameters {

    real a[N];

    real b[N];

    for (i in 1:N) {

        a[i] = a1 + a2 * G[i];

        b[i] = b1 + b2 * G[i];

    }

}

model {

    a1 ~ normal(0, 10);

    a2 ~ normal(0, 10);

    b1 ~ normal(0, 10);

    b2 ~ normal(0, 10);

    sgm ~ exponential(0.1);

    for (i in 1:N) {

        Y[i] ~ normal(a[i] + b[i] * X[i], sgm);

    }

}

 

表１のデータを分析するために、2つのExcelファイルを図１のように用意した。

図１

 

図１のファイルを読み込んで、リスト１のStanスクリプトにより分析するPythonスクリプトをリスト２のように用意した。温度Xを中心化した値Zを独立変数としていることに注意。

 

リスト２　ダミー変数を用いた回帰モデルの分析Pythonスクリプト（SweetTempDummy.py）

import pandas as pd

import numpy as np

import matplotlib.pyplot as plt

import pystan

import seaborn as sb

 

xlfile = pd.ExcelFile('fructose.xlsx')

print(xlfile)

df = pd.read_excel(xlfile)

print(df)

col_names = df.keys()

print(col_names)

Xf = df['Temperature'].values

Yf = df['Sweetness'].values

Gf = np.ones(len(Xf))         #    ダミー変数 G=1

print('Xf =\n', Xf)

print('Yf =\n', Yf)

print('Gf =\n', Gf)

 

xlfile = pd.ExcelFile('glucose.xlsx')

print(xlfile)

df = pd.read_excel(xlfile)

print(df)

col_names = df.keys()

print(col_names)

Xg = df['Temperature'].values

Yg = df['Sweetness'].values

Gg = np.zeros(len(Xf))        #   ダミー変数 G = 0

print('Xg =\n', Xg)

print('Yg =\n', Yg)

print('Gg =\n', Gg)

 

X = np.concatenate((Xf,Xg))

Y = np.concatenate((Yf, Yg))

G = np.concatenate((Gf, Gg))

print('X =\n', X)

print('Y =\n', Y)

print('G =\n', G)

 

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

print('Compilation ended.')

 

Z = X - np.mean(X)

Data = {'N':len(Y), 'Y':Y, 'X':Z, 'G':G}

fit = sm.sampling(data = Data, n_jobs = 1)

print(fit)

 

sgm_hat = np.median(fit['sgm'])

sb.kdeplot(fit['sgm'])

plt.title('Sigma: median = {0:.2}'.format(sgm_hat))

plt.savefig('FigSgm.png')

plt.show()

 

Y_hat = np.empty(len(Y))

a1_hat = np.median(fit['a1'])

a2_hat = np.median(fit['a2'])

b1_hat = np.median(fit['b1'])

b2_hat = np.median(fit['b2'])

 

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

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

plt.title('a1: median={0:.2f},   a2: median={1:.2f}'.format(

    a1_hat, a2_hat))

plt.legend()

plt.savefig('Figa1a2.png')

plt.show()

 

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

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

plt.title('b1: median={0:.4f},   b2: median={1:.4f}'.format(

    b1_hat, b2_hat))

plt.xticks(rotation = 45)

plt.legend()

plt.tight_layout()

plt.savefig('Figb1b2.png')

plt.show()

 

plt.plot(Z, Y, 'o')

print(np.min(Z), np.max(Z))

plt.plot([np.min(Z), np.max(Z)], [a1_hat + b1_hat*np.min(Z), a1_hat + b1_hat*np.max(Z)],

        label = 'glucose')

plt.plot([np.min(Z), np.max(Z)],

         [a1_hat + a2_hat + (b1_hat + b2_hat)*np.min(Z), a1_hat + a2_hat + (b1_hat + b2_hat)*np.max(Z)],

         label = 'fructose')

x_coord = np.arange(0, 61, 10)

z_coord = x_coord - np.mean(X)

plt.xticks(z_coord, x_coord)

plt.xlabel('$^{\circ}C$', fontsize = 14)

plt.ylabel('Sweetness')

plt.plot([0, 0], [np.min(Y), np.max(Y)], label = 'Mean of X = {0:.1f}'.format(np.mean(X)))

plt.legend()

plt.tight_layout()

plt.savefig('FigRegressDummy.png')

plt.show()

 

print('Mean of X = {0:.5f}'.format(np.mean(X)))

 

ck_b1 = [1 if v < 0 else 0 for v in fit['b1']]

print('Prop(b1 < 0) = {0:.5f}'.format(np.sum(ck_b1) / len(ck_b1)))

 

 

図１のデータファイルに対してリスト２を実行すると、以下の結果を得る。

 

まず、sigmaの事後分布が図２のように示される。

図２

 

続いて、a1とa2の事後分布が図３のように示される。ぽあらめーたa1およびa2ともに正の値であることが分かる。

図３

 

パラメータb1とb2の事後分布は図４のようである。

図４

 

パラメータb1の事後分布は０より左側に寄っているように見える。端末出力には

 

Prop(b1 < 0) = 0.91100

 

と出力されている。b1の事後分布の91％が負の側にあるので、０に近い値であるが、負の値であると言える。

 

データ点と回帰直線は図５のようである。赤の縦線は温度の平均値を示している。パラメータa1とa2の値は、この位置における温度と甘味度の関係を表している。

図５

 

グルコースの回帰直線が水平に近いことは、b1の値が０に近いことに対応している。b1のあたいは、０に近いが負の値であると考えられる事後分布であるので、温度とともにわずかであるが甘味度は減少していると言える。グルコースに比べてフラクトースの甘味度の変化は明瞭である。これはb2の事後分布（図４）がb1に比べて大きく負の側に移動していることに対応している。

 

 

引用文献

前橋健二（2011）．甘味の基礎知識．日本醸造協会誌、106、818-825．

 

 

Up