心理尺度の信頼性係数ω
ベイズ分析
心理尺度の信頼性係数としてオメガω(McDonald, 1999)は項目jの得点の因子モデル
に基づいて算出される(岡本、2014、第8章)。ここで、は測定対象の値(真値)、
は項目jに固有の特性など誤差項としてまとめられるものである。
は平均0、分散1に標準化されたものを考えるが、
は素点である。
項目が項目あるとし、それらの合計点
がテスト得点として算出されるとする。
である。このとき、オメガは次式で与えられる。
の平均が0、分散が1であることから式(3)は
と書ける。式番号(6.20a)はMcDonald (1999)におけるものである。
さらに、各が他の変数と独立であるとすると、式(6.20a)は
と書ける。
各項目の真値の係数(因子負荷量)
が同じであるとき、アルファ係数の公式が導かれる。
ベイズ分析例
ベイズ分析のためのStanスクリプトを用いる例をリスト1に、リスト1によるMCMCサンプルの分析を行うスクリプトをリスト2に示した。ファイル類は、まとめて圧縮ファイルomegapystan3files.zipとした。
Pystan3との組み合わせでArvizが使えるようになったので、更新版をこのウェブサイトに用意した(2025.03.22)。MCMCサンプルのRhat値とトレース図が得られるようにした。
リスト1のスクリプトは、WindowsにインストールされたUbuntu上で実行された。リスト2は、WindowsにインストールされたPythonで実行された。Windows上ではPyStan3は直接実行できないので、WindowsにWSL(Windows Subsystem for Linux)をインストールして、その上にUbuntuをインストールした。WindowsへのUbuntuのインストールはこのウェブサイトで説明している。PyStan3の簡単な説明のウェブサイトも用意している。WSLでは描画ができなかったので、リスト1で生成したMCMCサンプルの描画は、直接Windows上でPython 3.13により行った。リスト1のスクリプトファイルとリスト2のスクリプトファイルを同じフォルダに置いておくと、リスト1のUbuntuにおける実行に続いてリスト2をWindowsで実行することがスムーズに行うことができる。
リスト1のStanスクリプトは、式(1)
のモデルに対するものであるが、式(1)におけるTはリスト1ではFで表されている。Tはtrue valueを意味しているが、Fはfactor scoreであり、すなわち真値である。
パラメータFは、平均0,標準偏差1の標準正規分布に従うという仮定が置かれているが、事前分布での制約ではshrinkageが起こり、このためが過大評価されるという傾向が見られる。これを避けるために、Fの標準化を行っている。
また、MCMCサンプリングにおいてlabel
switchingを避けるためにパラメータに非負の制約を置いている。
事前分布は、soft constraintを用いたが、Gelmanら(2021、p. 124)より緩いものを用いた。これは、Gelmanらでは、データ数がsparseの場合にも対応できるように設定されているが、尺度の信頼性を求めるときのデータはある程度の多さがあると考えられるからである。データ数が少なくてパラメータ推定が不安定な場合は、それに応じて事前分布を調整すればよい。PyStan 3では、事前分布として分散の大きいものを設定すると警告が出るので、一様分布に近いものとして分散の大きい分布を用いるときは、まず警告の出ない事前分布のパラメータを用意して、これを変換して分散の大きい分布としている(rescaling)。
リスト1はPyStan 3対応である。PyStan 3をWindowsで利用するときは、WindowsにUbuntuをインストールしてUbuntu上で利用するが、UbuntuのWindowsへのインストール、およびUbuntuでのPythonの準備は簡単なコマンドで出来る。このウェブサイトで説明しているので、必要ならば参照されたい。
リスト1のスクリプトを実行すると、入力データファイル名を聞いてくる。
(py312)
*****/omegapystan3files$ python EstOmega_ps3.py
Input Excel
data file (*.xlsx) =
入力データファイルはExcelファイルとして図1のように用意する。なお、図1のデータの再分析もこのウェブサイトの後半に示している。
図1
1行目に変数名を設定して、次の2行目からデータを設定する。ただし、1列目はデータの識別用値を設定する。1列目の値は、分析では無視される。
図1のファイルを入力データするときは、入力データファイル名を次のように設定する。
Input Excel
data file (*.xlsx) = Data10PS3.xlsx
ファイル名を設定すると、計算が始まる。MCMCサンプリングが終了すると、出力ファイルResults.txtが保存されたことを表示してスクリプトの実行が終了する。
Mu[1]: Q1 = 50.773, Med = 51.782, Q3 = 52.784
Mu[2]: Q1 = 53.389, Med = 54.234, Q3 = 55.089
Mu[3]: Q1 = 49.044, Med = 50.187, Q3 = 51.332
Mu[4]: Q1 = 51.570, Med = 52.644, Q3 = 53.673
Mu[5]: Q1 = 51.922, Med = 52.893, Q3 = 53.820
Mu[6]: Q1 = 44.403, Med = 46.061, Q3 = 47.776
Mu[7]: Q1 = 51.574, Med = 52.507, Q3 = 53.406
Mu[8]: Q1 = 51.353, Med = 52.288, Q3 = 53.155
Mu[9]: Q1 = 42.242, Med = 43.764, Q3 = 45.282
Mu[10]: Q1 = 53.064, Med = 54.047, Q3 = 54.995
Results.txt was saved.
(py312) *****/omegapystan3files$
出力ファイルResults.txtは、エディタで開くことができる。
リスト1のスクリプトの実行中に保存された他のファイルは、次のリスト2のスクリプトの実行で読み込まれて分析される。リスト1とリスト2は同じフォルダ内に置かれているものとする。
リスト2のスクリプトを実行すると、まずの事後分布のグラフが表示される(図2)。
図2
図2のWindowを閉じると、の事後分布が表示される(図3)。
図3
図3のWindowを閉じると、の事後分布のグラフが表示される(図4)。
図4
図4のWindowを閉じると、次に、変数ごとに各変数における真値の分散の各変数の分散に対する比
の事後分布のグラフが表示される(図5)。この値は、各変数に占める真値の割合と考えられる。
図5
図5のWindowを閉じると、ωの事後分布のグラフが表示される(図6)。
図6
グラフのタイトル部に、事後分布の4分位数が表示されている。点推定値としては、MAP推定値(モード)がよいという報告(岡本、2022)があるが、計算時間が掛かるので、ここではGelmanら(2021)に従って、中央値を点推定値とした。
図6のフォームを閉じると、図7の散布図が表示される。
図7
横軸に各の推定値
(事後分布の中央値)、縦軸にテストスコアの和
をとり、プロットしたものである。点
がほぼ直線沿いに散らばっていることが判る。
図7のWindowを閉じると、スクリプトの実行終了となる。端末画面には、次のメッセージが表示されている。
Results_1.txt
was saved.
出力ファイルResults_1.txtを開くと、以下のようになっている。
data...
[[ 1 79 46 68 65 89 38 58 70 49 65]
[ 2 48 71 61 50 43 27 48 54 34 61]
[ 3 38 42 38 45 55 43 47 39 65 54]
・
・
・
[48 36 46 49 54 31 75 43 40 30 49]
[49 44 24 56 51 58 34 42 52 67 48]
[50 71 75 64 84 58 25 67 65 47 56]]
Mu...
Mu[1]: Q1 = 50.830, Med = 51.831, Q3 = 52.865
Mu[2]: Q1 = 53.353, Med = 54.220, Q3 = 55.137
Mu[3]: Q1 = 49.028, Med = 50.159, Q3 = 51.319
Mu[4]: Q1 = 51.609, Med = 52.646, Q3 = 53.687
Mu[5]: Q1 = 51.915, Med = 52.871, Q3 = 53.848
Mu[6]: Q1 = 44.512, Med = 46.133, Q3 = 47.731
Mu[7]: Q1 = 51.613, Med = 52.482, Q3 = 53.394
Mu[8]: Q1 = 51.407, Med = 52.248, Q3 = 53.139
Mu[9]: Q1 = 42.293, Med = 43.748, Q3 = 45.241
Mu[10]: Q1 = 53.102, Med = 54.068, Q3 = 55.051
Sigma...
Sigma[1]: Q1 = 9.642 Med = 10.410, Q3 = 11.247
Sigma[2]: Q1 = 8.284 Med = 8.978, Q3 = 9.735
Sigma[3]: Q1 = 10.983 Med = 11.842, Q3 = 12.747
Sigma[4]: Q1 = 9.919 Med = 10.687, Q3 = 11.617
Sigma[5]: Q1 = 9.202 Med = 9.970, Q3 = 10.813
Sigma[6]: Q1 = 16.324 Med = 17.401, Q3 = 18.750
Sigma[7]: Q1 = 8.567 Med = 9.282, Q3 = 10.048
Sigma[8]: Q1 = 8.575 Med = 9.357, Q3 = 10.224
Sigma[9]: Q1 = 14.564 Med = 15.605, Q3 = 16.794
Sigma[10]: Q1 = 9.443 Med = 10.207, Q3 = 11.154
Lambda...
Lambda[1]: Q1 = 12.019 Med = 13.024, Q3 = 14.026 (Item_1)
Lambda[2]: Q1 = 11.692 Med = 12.602, Q3 = 13.466 (Item_2)
Lambda[3]: Q1 = 11.536 Med = 12.743, Q3 = 13.908 (Item_3)
Lambda[4]: Q1 = 15.121 Med = 16.120, Q3 = 17.180 (Item_4)
Lambda[5]: Q1 = 12.310 Med = 13.275, Q3 = 14.242 (Item_5)
Lambda[6]: Q1 = 0.537 Med = 1.211, Q3 = 2.216 (Item_6)
Lambda[7]: Q1 = 11.363 Med = 12.291, Q3 = 13.136 (Item_7)
Lambda[8]: Q1 = 14.754 Med = 15.739, Q3 = 16.642 (Item_8)
Lambda[9]: Q1 = 1.191 Med = 2.267, Q3 = 3.587 (Item_9)
Lambda[10]: Q1 = 14.841 Med = 15.834, Q3 = 16.835 (Item_10)
Lambda_j**/Var(X_j)...
Lambda[1]**2 / Var(X[1]): Q1 = 0.555, Med = 0.610, Q3 = 0.662
Lambda[2]**2 / Var(X[2]): Q1 = 0.610, Med = 0.662, Q3 = 0.710
Lambda[3]**2 / Var(X[3]): Q1 = 0.474, Med = 0.538, Q3 = 0.594
Lambda[4]**2 / Var(X[4]): Q1 = 0.644, Med = 0.695, Q3 = 0.738
Lambda[5]**2 / Var(X[5]): Q1 = 0.584, Med = 0.640, Q3 = 0.688
Lambda[6]**2 / Var(X[6]): Q1 = 0.001, Med = 0.005, Q3 = 0.016
Lambda[7]**2 / Var(X[7]): Q1 = 0.580, Med = 0.636, Q3 = 0.688
Lambda[8]**2 / Var(X[8]): Q1 = 0.691, Med = 0.739, Q3 = 0.778
Lambda[9]**2 / Var(X[9]): Q1 = 0.006, Med = 0.021, Q3 = 0.050
Lambda[10]**2 / Var(X[10]): Q1 = 0.657, Med = 0.705, Q3 = 0.748
omega...
Q1 = 0.898, Med = 0.905, Q3 = 0.911
まず、入力データが出力されている。
続いて、MCMCサンプリングの統計量が出力されている。
事後分布などのグラフで表示されていた情報が、統計量の数値データとして読むことができる。
図4()、図5(
)における他とは著しく異なる項目は、上の統計量から6番目と9番目の項目であることが解る。6番目と9番目を除いて再分析を行う。
図1のデータから、6番目と9番目の項目のデータを抜いたものを、図A1に示すように用意した。
図A1 図1のデータから6番目と9番目のデータを抜いたもの
図A1のデータに対して、まずリスト1のスクリプトでMCMCサンプリングを行った。
(py312) *****/omegapystan3files$
python EstOmega_ps3.py
Input Excel data file (*.xlsx) =
Data10PS3_69.xlsx
リスト1のスクリプトによるMCMCサンプリングの後、MCMCサンプルの分析をリスト2のスクリプトで行った。
の事後分布は、図A2のようである。
図A2
の事後分布は図A3のようである。
図A3
の事後分布は、図A4のようである。
図A4
の事後分布は図A5のようである。
図A5
ωの事後分布は図A6のようである。図6のものより右に寄っていることが分かる。
図A6
点の散布図は図A7のようである。図7のものより、点が狭い帯内にまとまっていることが分かる。
図A7
出力ファイルResults_1.txtを開くと、以下のようになっている。
data...
[[ 1 79 46 68 65 89 58 70 65]
[ 2 48 71 61 50 43 48 54 61]
[ 3 38 42 38 45 55 47 39 54]
・
・
・
[48 36 46 49 54 31 43 40 49]
[49 44 24 56 51 58 42 52 48]
[50 71 75 64 84 58 67 65 56]]
Mu...
Mu[1]: Q1 = 50.764, Med = 51.760, Q3 = 52.738
Mu[2]: Q1 = 53.410, Med = 54.258, Q3 = 55.118
Mu[3]: Q1 = 49.034, Med = 50.143, Q3 = 51.341
Mu[4]: Q1 = 51.568, Med = 52.626, Q3 = 53.648
Mu[5]: Q1 = 51.974, Med = 52.900, Q3 = 53.838
Mu[6]: Q1 = 51.618, Med = 52.488, Q3 = 53.403
Mu[7]: Q1 = 51.372, Med = 52.248, Q3 = 53.125
Mu[8]: Q1 = 53.082, Med = 54.033, Q3 = 54.994
Sigma...
Sigma[1]: Q1 = 9.636 Med = 10.431, Q3 = 11.266
Sigma[2]: Q1 = 8.233 Med = 8.907, Q3 = 9.690
Sigma[3]: Q1 = 11.030 Med = 11.894, Q3 = 12.825
Sigma[4]: Q1 = 9.887 Med = 10.713, Q3 = 11.659
Sigma[5]: Q1 = 9.228 Med = 9.968, Q3 = 10.830
Sigma[6]: Q1 = 8.638 Med = 9.297, Q3 = 10.042
Sigma[7]: Q1 = 8.619 Med = 9.381, Q3 = 10.245
Sigma[8]: Q1 = 9.387 Med = 10.177, Q3 = 11.057
Lambda...
Lambda[1]: Q1 = 11.991 Med = 13.022, Q3 = 14.071 (Item_1)
Lambda[2]: Q1 = 11.698 Med = 12.601, Q3 = 13.507 (Item_2)
Lambda[3]: Q1 = 11.491 Med = 12.644, Q3 = 13.867 (Item_3)
Lambda[4]: Q1 = 15.017 Med = 16.109, Q3 = 17.154 (Item_4)
Lambda[5]: Q1 = 12.293 Med = 13.268, Q3 = 14.273 (Item_5)
Lambda[6]: Q1 = 11.346 Med = 12.275, Q3 = 13.163 (Item_7)
Lambda[7]: Q1 = 14.737 Med = 15.733, Q3 = 16.612 (Item_8)
Lambda[8]: Q1 = 14.813 Med = 15.798, Q3 = 16.847 (Item_10)
Lambda_j**/Var(X_j)...
Lambda[1]**2 / Var(X[1]): Q1 = 0.551, Med = 0.610, Q3 = 0.663
Lambda[2]**2 / Var(X[2]): Q1 = 0.610, Med = 0.667, Q3 = 0.715
Lambda[3]**2 / Var(X[3]): Q1 = 0.468, Med = 0.530, Q3 = 0.589
Lambda[4]**2 / Var(X[4]): Q1 = 0.640, Med = 0.693, Q3 = 0.738
Lambda[5]**2 / Var(X[5]): Q1 = 0.584, Med = 0.640, Q3 = 0.690
Lambda[6]**2 / Var(X[6]): Q1 = 0.577, Med = 0.637, Q3 = 0.685
Lambda[7]**2 / Var(X[7]): Q1 = 0.689, Med = 0.737, Q3 = 0.778
Lambda[8]**2 / Var(X[8]): Q1 = 0.658, Med = 0.706, Q3 = 0.750
omega...
Q1 = 0.931, Med =
0.936, Q3 = 0.941
他と特に互いに異なる値の項目は認められない。ωの値も高くなっており
10項目のときの
Q1 = 0.898, Med = 0.905, Q3 = 0.911
から、第6項目と第9項目を抜いた8項目でのωは
Q1 = 0.931, Med =
0.936, Q3 = 0.941
であり、四分位数の比較でも8項目に減らしたときのωの方が高い値であると言えよう。
import numpy as np
import pandas as pd
import stan
import pickle
StanCode = """
//
//
Yasuharu Okamoto, 2023.08, 2025.03
//
data {
int
N;
int
M;
array[N,M] real X;
vector[M] i_mu;
vector[M] i_std;
}
parameters {
vector[N] vF;
vector[M] v_mu;
vector<lower=0.0>[M] v_Lmbd;
vector<lower=0.0>[M] v_sgm;
}
transformed parameters {
vector[N] F;
real
sumF;
real
ssF;
real
meanF;
real
sdF;
vector[M] mu;
vector[M] Lmbd;
vector[M] sgm;
array[N] vector[M] V;
//
//
Standardized F values
//
sumF
= 0.0;
ssF
= 0.0;
for
(i in 1:N) {
sumF += vF[i];
ssF += pow(vF[i], 2);
}
meanF =
sumF / N;
sdF
= pow((ssF - N * pow(meanF, 2)) / (N - 1), 0.5);
F =
(vF - meanF) / sdF; // 標準化
//
//
Rescaled parameters
//
mu =
i_mu + 10.0 .* i_std .* v_mu;
Lmbd
= 10.0 .* i_std .* v_Lmbd;
sgm
= 0.001 + 10.0 .* i_std .* v_sgm;
//
//
Regression model
//
for
(i in 1:N)
V[i] = mu + Lmbd * F[i];
}
model {
vF ~
normal(0.0, 1.0);
v_mu
~ normal(0.0, 10.0);
v_Lmbd ~ exponential(0.1);
v_sgm ~ exponential(0.1);
for
(i in 1:N)
X[i] ~ normal(V[i], sgm);
}
generated quantities {
real
omg;
// omega
vector[M] L2V;
real
v1;
real
v2;
v1 =
0.0;
v2 =
0.0;
L2V
= pow(Lmbd, 2.0) ./ (pow(Lmbd, 2.0) + pow(sgm, 2.0));
for
(j in 1:M) {
v1 += Lmbd[j];
v2 += pow(sgm[j], 2);
}
omg
= pow(v1, 2) / (pow(v1, 2) + v2);
}
"""
pd.set_option('display.max_columns',
10000)
pd.set_option('display.max_rows',
10000)
fout = open('Results.txt', 'w')
#
# Excelデータファイルの読み込み
#
s = input("Input Excel data file
(*.xlsx) = ")
xlsx_fl = pd.ExcelFile(s)
data_xlsx = pd.read_excel(xlsx_fl)
with open('data.pkl', 'wb') as f:
pickle.dump(data_xlsx, f)
var_names = [v for v in
data_xlsx.keys()][1:]
print('var_names =', var_names)
data = data_xlsx.values
print('data...\n', data)
fout.write(f'Input data file = {s}\n')
fout.write(f'\ndata...\n{data}\n')
X = data[:,1:]
mus = np.mean(X, axis=0)
stds = np.std(X, axis=0)
Data = {'N': len(X), 'M': len(X[0]),
'X': X, 'i_mu': mus, 'i_std': stds}
sm = stan.build(StanCode, data = Data)
fit = sm.sample()
df = fit.to_frame() # DataFrame型への変換
with open('fit.pkl', 'wb') as f:
pickle.dump(df, f)
print()
fout.write('\n\nMu...\n')
for j in range(len(X[0])):
q1, med, q3 =
np.percentile(df[f'mu.{j+1}'], [25, 50, 75])
print(f'Mu[{j+1}]:', end =
'')
print(f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}')
fout.write(f'Mu[{j+1}]: Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}\n')
fout.close()
print('\nResults.txt was saved.\n')
import numpy as np
import pandas as pd
import pickle
import matplotlib.pyplot as plt
import seaborn as sb
fout = open('Results_1.txt', 'w')
with open('data.pkl', 'rb') as f: # 保存されたデータの読み込み
data_xlsx = pickle.load(f)
var_names = [v for v in
data_xlsx.keys()][1:]
print('var_names =', var_names)
data = data_xlsx.values
print('data...\n', data)
fout.write(f'\ndata...\n{data}\n')
X = data[:,1:]
mus = np.mean(X, axis=0)
stds = np.std(X, axis=0)
with open('fit.pkl', 'rb') as f: # 保存されたMCMCサンプルの読み込み
df = pickle.load(f)
#
# パラメータmuの4分位数の計算
#
print()
fout.write('\n\nMu...\n')
for j in range(len(X[0])):
q1, med, q3 =
np.percentile(df[f'mu.{j+1}'], [25, 50, 75])
print(f'Mu[{j+1}]:', end =
'')
print(f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}')
fout.write(f'Mu[{j+1}]: Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}\n')
#
# muの事後分布
#
for j in range(len(X[0])):
sb.kdeplot(df[f'mu.{j+1}'],
label=f'{var_names[j]}')
plt.title(r'$\mu$', fontsize = 16)
plt.xlabel('')
plt.legend(bbox_to_anchor=(1.0, 1.0))
plt.tight_layout()
plt.show()
#
#
Sigmaの4分位数の計算
#
print()
fout.write('\nSigma...\n')
for j in range(len(X[0])):
q1, med, q3 =
np.percentile(df[f'sgm.{j+1}'], [25, 50, 75])
print(fr'Sigma[{j+1}]:', end
= '')
print(f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}')
fout.write(f'Sigma[{j+1}]:
Q1 = {q1:.3f}' +
f' Med = {med:.3f}, Q3 = {q3:.3f}\n')
#
#
Sigmaの事後分布
#
for j in range(len(X[0])):
sb.kdeplot(df[f'sgm.{j+1}'],
label=f'{var_names[j]}')
plt.title(r'$\sigma$', fontsize = 16)
plt.xlabel('')
plt.legend(bbox_to_anchor=(1.0, 1.0))
plt.tight_layout()
plt.show()
#
# Lambdaの四分位数の計算
#
print()
fout.write('\nLambda...\n')
for j in range(len(X[0])):
q1, med, q3 =
np.percentile(df[f'Lmbd.{j+1}'], [25, 50, 75])
print(fr'Lambda[{j+1}]:',
end = '')
print(f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}', end = '')
print(f' ({var_names[j]})')
fout.write(f'Lambda[{j+1}]:
Q1 = {q1:.3f}' +
f' Med = {med:.3f}, Q3 = {q3:.3f}' +
f' ({var_names[j]})\n')
#
# Lambdaの事後分布
#
for j in range(len(X[0])):
sb.kdeplot(df[f'Lmbd.{j+1}'], label=f'{var_names[j]}')
plt.title(r'$\lambda$', fontsize = 16)
plt.xlabel('')
plt.legend(bbox_to_anchor=(1.0, 1.0))
plt.tight_layout()
plt.savefig('FigLmbd.png')
plt.show()
#
#
真値の割合の計算
#
print()
fout.write('\nLambda_j**/Var(X_j)...\n')
for j in range(len(X[0])):
q1, med, q3 =
np.percentile(df[f'L2V.{j+1}'], [25, 50, 75])
print(f'Lambda[{j+1}]**2 /
Var(X[{j+1}]):', end = '')
print(f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}')
fout.write(f'Lambda[{j+1}]**2 / Var(X[{j+1}]):' +
f' Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}\n')
for j in range(len(X[0])):
sb.kdeplot(df[f'L2V.{j+1}'],
label=f'{var_names[j]}')
plt.title(r'$\lambda^{2}_{j}/Var(X_{j})$',
fontsize = 16)
plt.xlabel('')
plt.legend(bbox_to_anchor=(1.0, 1.0))
plt.tight_layout()
plt.show()
#
#
Omegaの四分位数の計算
#
print()
q1, med, q3 =
np.percentile(df[f'omg'], [25, 50, 75])
print(f'omega: Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}')
fout.write('\nomega...\n')
fout.write(f'
Q1 = {q1:.3f}, Med =
{med:.3f}, Q3 = {q3:.3f}\n')
#
#
Omegaの事後分布
#
sb.kdeplot(df['omg'])
plt.title(f'Q1 = {q1:.3f}, Med = {med:.3f}, Q3 = {q3:.3f}', fontsize = 14)
plt.xlabel(r'$\omega$', fontsize = 16)
plt.tight_layout()
plt.show()
#
# (因子得点F、尺度得点和)の散布図
#
sumXs = np.sum(X, axis = 1)
Fs = []
for i in range(len(X)):
Fs.append(np.median(df[f'F.{i+1}']))
plt.plot(Fs, sumXs, 'o')
plt.xlabel(r'$\widehat{F}_{median}$',
fontsize = 14)
plt.ylabel(r'$\sum{X_{j}}$', fontsize
= 14)
plt.tight_layout()
plt.show()
fout.close()
print('\nResults_1.txt was saved.\n')
Gelman, A., Hill, J., & Vehtari, A. (2021). Regression and Other Stories. Cambridge University Press.
McDonald, R. P. (1999). Test Theory: A unified treatment. Psychology Press.