心理尺度の信頼性係数ω
尺度の分析を、主因子法のアルゴリズムを適用する場合、直接最小2乗法を適用して特異値分解を利用する場合、ベイズ分析を用いる場合を扱う。
主因子法は途中の繰り返し計算が収束しないことがあるが、特異値分解は常に解が得られるようである。
ベイズ分析では正規分布が仮定されているが、特異値分解による方法の場合は平均が0、分散が1という位置とサイズについての仮定だけで、分布の形状についての仮定はない。
心理尺度の信頼性係数としてオメガω(McDonald,
1999)は項目jの得点の因子モデル
に基づいて算出される(岡本、2014、第8章)。ここで、は測定対象の値(真値)、
は項目jに固有の特性など誤差項としてまとめられるものである。
は平均0、分散1に標準化されたものを考えるが、
は素点である。
項目が項目あるとし、それらの合計点
がテスト得点として算出されるとする。
である。このとき、オメガは次式で与えられる。
の平均が0、分散が1であることから式(3)は
と書ける。式番号(6.20a)はMcDonald (1999)におけるものである。
さらに、各が他の変数と独立であるとすると、式(6.20a)は
と書ける。
各項目の真値の係数(因子負荷量)
が同じであるとき、アルファ係数の公式が導かれる。
モデル(1)におけるの推定を、主因子法のアルゴリズムで行う。すなわち、因子分析における主因子法の計算方法を共分散行列に適用して因子負荷量
を推定しする。詳しくは、<岡本安晴「心理学データ分析と測定」勁草書房>の「第8章 心理学的特性を測る:心理尺度」を参照されたい。ベイズ分析のスクリプト例は、このウェブサイトの後半に示した。
主因子法のアルゴリズムを用いるコードは、以下のように簡単である。
psi2 =
numpy.zeros(n_vars)
A =
numpy.zeros(n_vars)
while True:
A0 = A.copy()
C = S.copy()
for i in range(n_vars):
C[i][i] -= psi2[i]
w, vctr = numpy.linalg.eigh(C)
A = vctr[:, -1] *
(w[-1] ** 0.5)
for i in range(n_vars):
psi2[i] = S[i][i] - (A[i] ** 2)
diff = ((A0 - A) ** 2).sum()
if diff < 1.0e-7:
break
プログラムはサンプルデータファイルとともに圧縮ファイルfilesEstOmega.zipにまとめた。ダウンロード・展開すれば自由に利用できる。
なお、Jupyter Notebook用のプログラムも用意した。圧縮ファイルfilesEstOmega_jn.zipにまとめている。Jupyter Notebookでの実行画面例をここに示した。
プログラムはPythonで書かれているが、Pythonの簡単なインストール法と使用法を説明したウェブサイトも用意している。
圧縮ファイルをダウンロード展開して作成されるフォルダにあるプログラムEstOmegaMain.pyを実行すると、入力データファイル名の入力が求められる(図1)。
図1
入力データファイルはExcelファイルとして用意する。図2は、10項目Item_1からItem_10に対する50人分のデータを表している。
図2.a 入力データファイル例(Data10.xlsx)
図2.b 入力データファイル例(Data10.xlsx)
1行目には変数名が書かれている。第1列目はデータの識別値であり、数値でも文字列でも何でもよい。識別値であるが、同じ値が複数回用いられてもよい。
2行目には、その変数が信頼性係数の計算において含めるかどうかを示すものである。整数値1が設定されている変数が尺度を構成する変数として含められ、整数値0が設定されていると尺度から除かれる。一部の変数を対象とする分析例は、ウェブサイトの後半で説明している。
データ値は3行目以降に書かれる。データ値は、逆転項目が用いられている場合は、逆転項目の値を逆転して他の項目と値の方向を合わせておく必要がある。
入力データファイル名を設定してEnterキーを押すと、次に出力用テキストデータファイル名の設定が求められる(図3)。
図3
図3では、入力データファイル名として図2のファイルData10.xlsxが設定され、出力用テキストファイル名としてResults.txtが設定されている。出力用テキストファイル名は任意であるが、ファイル拡張子が「.txt」であることに注意。
出力用テキストファイル名を入力すると計算が始まる。
計算が終了すると、出力用テキストファイルが保存されたことを示すメッセージが表示されて、プログラムの実行終了となる。
図4
実行終了後、出力用テキストファイルを開くことができる。
出力ファイルをPython IDLEで開こうとすると、図5のダイアログボックスが表示されることがある。
図5
この場合は、「s-jis」を設定して(図6)、「OK」ボタンをクリックすればよい。
図6
出力用テキストファイルの内容は、以下のようになっている。説明は、このリストの後で行っている。
Input data file = Data10.xlsx
n_data =50
n_vars = 10
var_names: ['Item_1', 'Item_2',
'Item_3', 'Item_4', 'Item_5', 'Item_6', 'Item_7', 'Item_8', 'Item_9',
'Item_10']
α = 0.879
Lambda[1]**2 / Var(X[1]) = 0.620
Lambda[2]**2 / Var(X[2]) = 0.682
Lambda[3]**2 / Var(X[3]) = 0.553
Lambda[4]**2 / Var(X[4]) = 0.702
Lambda[5]**2 / Var(X[5]) = 0.664
Lambda[6]**2 / Var(X[6]) = 0.010
Lambda[7]**2 / Var(X[7]) = 0.652
Lambda[8]**2 / Var(X[8]) = 0.746
Lambda[9]**2 / Var(X[9]) = 0.011
Lambda[10]**2 / Var(X[10]) = 0.714
Item_1: Lambda[1] = 12.990
Item_2: Lambda[2] = 12.642
Item_3: Lambda[3] = 12.765
Item_4: Lambda[4] = 16.076
Item_5: Lambda[5] = 13.397
Item_6: Lambda[6] = -1.719
Item_7: Lambda[7] = 12.301
Item_8: Lambda[8] = 15.702
Item_9: Lambda[9] = 1.611
Item_10: Lambda[10] = 15.803
GFI = 0.995
ω = 0.907
入力データファイル名、データ数、変数名に続いて、信頼性係数アルファαが出力されている。このαの値は、McDonald (1999)がGuttman-Cronbach alphaとも呼んでいるもので、いわゆるアルファ係数である。アルファ係数については、<岡本安晴「心理学データ分析と測定」勁草書房>のpp. 126-128で説明している。
続いて、変数ごとに各変数における真値の分散の各変数の分散に対する比
の値が出力されている。これは、各変数に占める真値の割合と考えられる。
上の出力では、6番目と9番目の値が低い。これらは尺度から除くのが望ましい変数と考えられる。
次に、各の値が出力されている。モデル(1)が素点
に対する因子モデルであるので、因子負荷量
の値は1より大きくなり得る。6番目の項目の
の値が負であるが、これは変数の値の動きが真値と逆であることを意味しているので、逆転項目の可能性がある。また、値の絶対値も小さいので尺度から除くことも考えられる。
次に、GFIが示されているが、これはモデル(1)の適合度(goodness of fit index)を表し(McDonald (1999)、GFIは1に近いほど適合度が高く、0に近いほど適合度が低い(岡本安晴「心理学データ分析と測定」勁草書房の「第8章 心理学的特性を測る」、pp. 124-125参照)。
最後にの値が出力されている。式(6.20a)による値である。
信頼性係数アルファが
であるのに対し、信頼性係数ωは
である。
信頼性係数αが信頼性係数ωより小さい値であることが確認できる。
上の出力におけるの値を見ると、項目Item-6の値は負であり、Item-9の値は他と比べて小さい。これらの項目を除いて以下のように再計算を行った。
Item-6とItem-9をアルファ係数の算出から除くため、上の入力用ファイルを図A.1に示すように変更する。
図A.1 項目Item_6とItem_9を除外(ファイル名Data10_69.xlsx)
項目Item_6と項目Item_9の2行目のセルに0を設定することにより、データ識別変数DataIDの値とともに分析対象の変数から除外する。
プログラム実行時の入力データファイルとして、図A.1のファイルData10_69.xlsxを指定する(図A.2)。
図A.2
このときの出力ファイルは、以下のようになる。
Input data file = Data10_69.xlsx
n_data =50
n_vars = 8
var_names: ['Item_1', 'Item_2',
'Item_3', 'Item_4', 'Item_5', 'Item_7', 'Item_8', 'Item_10']
α = 0.940
Lambda[1]**2 / Var(X[1]) = 0.618
Lambda[2]**2 / Var(X[2]) = 0.683
Lambda[3]**2 / Var(X[3]) = 0.546
Lambda[4]**2 / Var(X[4]) = 0.709
Lambda[5]**2 / Var(X[5]) = 0.659
Lambda[6]**2 / Var(X[6]) = 0.653
Lambda[7]**2 / Var(X[7]) = 0.749
Lambda[8]**2 / Var(X[8]) = 0.715
Item_1: Lambda[1] = 12.968
Item_2: Lambda[2] = 12.648
Item_3: Lambda[3] = 12.685
Item_4: Lambda[4] = 16.145
Item_5: Lambda[5] = 13.350
Item_7: Lambda[6] = 12.306
Item_8: Lambda[7] = 15.736
Item_10: Lambda[8] = 15.817
GFI = 0.999
ω = 0.942
信頼性係数ωは
となり、すべての変数を対象とした最初の信頼積係数ω
より、大きくなっている。
また、各変数における真値の分散の変数の分散に対する比
Lambda[j]**2 / Var(X[j])
に特に小さい値も見当たらない。
主因子法のアルゴリズムの場合、繰り返し計算が収束しない場合がある。これに対して、特異値分解は常に成功することが期待される。
において、特異値分解に定数項の影響が負荷量
の推定に影響を与えるのを避けるため、
いま、個人の得点を
で、特性値
を
で表し、モデル(PC.1)を次式で表す。Pythonスクリプトでは、Tは紛らわしいので避けた次第である。
モデル(PC.2)において、誤差項を最小2乗の意味で最小にする
と
を求める。
いま、の特異値分解を次式の形で表す。
また、
と置けば、
より、
を得る。したがって、
と置けば、は平均0,分散
であり、
と
は式(PC.2)の最小2乗解である。因子
の平均が0,分散が1ということは、尺度の原点と単位を決めるもので、本質的な制約ではない。
上のモデルに基づくプログラムをリストSVD1のように用意した。プログラムファイルなどは圧縮ファイルfilesEstOmegaSVD.zipとしてまとめた。
なお、プログラムはJupyter Notebook用も作成して、圧縮ファイルfilesEstOmegaSVD_jn.zipとしてまとめた。Jupyter Notebookでの実行例の画面はここに用意した。
リストSVD1を実行すると、入力データファイル名の入力が求められる。
Input Excel
data file (*.xlsx) =
入力データファイルは、Excelファイル(*.xlsx)として図SVD1のように用意する。
図SVD1
1行目に変数名を置く。2行目は変数の選択を示す。1であれば分析に含められるが、1以外の変数は分析から除外される。図SVD1では、変数DataID、Item_6、Item_9が分析から除かれている。
3行目以降にデータ値を書く。
次の例では、図SVD1のExcelファイル名が設定されている。入力データファイル名を入力すると、次に出力ファイル名の設定が求められる。出力ファイル名は任意であるが、テキストファイルなので、ファイル名拡張子を「.txt」とする。
Input Excel data file (*.xlsx) =
Data10_69.xlsx
Output text data file (*.txt) =
次の例では、Results.txtと入力している。
Input Excel data file (*.xlsx) =
Data10_69.xlsx
Output text data file (*.txt) =
Results.txt
出力ファイル名を入力すると、計算が始まり、特異値分解における特異値のスクリープロットが表示される(図SVD2)。
図SVD2
1番目の特異値が特に大きい値であり、2番目以降はなだらかに減少している。モデル(1)、(2)の1因子モデルが成り立っていると考えられる。
図SVD2を閉じると、次に図SVD3の散布図が表示される。
図SVD3
モデルによる予測値
と、データ値
の散布図である。
図SVD3を閉じると、プログラムの終了である。
端末には、以下のようなメッセージが表示されている。
Results_Scree.png was saved.
ResultsScatter.png was saved.
Results.txt was saved.
スクリープロットがResults_Scree.pngに、散布図がResultsScatter.pn保存されていることを意味している。これらの名前は、出力ファイル名を基に作成されている。
出力ファイルは、以下のような内容である。
Input data file = Data10_69.xlsx
n_data =50
n_vars = 8
var_names: ['Item_1', 'Item_2',
'Item_3', 'Item_4', 'Item_5', 'Item_7', 'Item_8', 'Item_10']
Item_1 ==> X[1]
Item_2 ==> X[2]
Item_3 ==> X[3]
Item_4 ==> X[4]
Item_5 ==> X[5]
Item_7 ==> X[6]
Item_8 ==> X[7]
Item_10 ==> X[8]
Lambda[1]**2 / Var(X[1]) = 0.648
Lambda[2]**2 / Var(X[2]) = 0.690
Lambda[3]**2 / Var(X[3]) = 0.595
Lambda[4]**2 / Var(X[4]) = 0.747
Lambda[5]**2 / Var(X[5]) = 0.682
Lambda[6]**2 / Var(X[6]) = 0.665
Lambda[7]**2 / Var(X[7]) = 0.766
Lambda[8]**2 / Var(X[8]) = 0.746
Item_1: Lambda[1] = 13.280
Item_2: Lambda[2] = 12.717
Item_3: Lambda[3] = 13.244
Item_4: Lambda[4] = 16.572
Item_5: Lambda[5] = 13.578
Item_7: Lambda[6] = 12.422
Item_8: Lambda[7] = 15.906
Item_10: Lambda[8] = 16.154
GFI = 0.997
ω = 0.999
分析された変数名が、プログラムでの変数として何番目であるかの一覧が示されている。
最後に、ωの値も書き出されている。
import numpy
import pandas as pd
import matplotlib.pyplot as plt
s = input("Input Excel data file
(*.xlsx) = ")
xlsx_fl = pd.ExcelFile(s)
data_xlsx =
pd.read_excel(xlsx_fl, header =
None)
data = data_xlsx.values
s_out = input("Output text data
file (*.txt) = ")
fout = open(s_out, "w")
fout.write("Input data file =
" + s + "\n")
data = data.T
n_vars = 0
var_names = []
X = []
for v in data:
if v[1] == 1:
var_names.append(v[0])
X.append(list(v[2:]))
n_vars += 1
n_data = len(X[0])
print('n_data =', n_data)
print('n_vars =', n_vars)
print('var_names = :', var_names)
for i, nm in enumerate( var_names):
print(f'{nm:<.10s} ==>
X[{i+1}]')
print()
fout.write('\nn_data
={}\n'.format(n_data))
fout.write('n_vars =
{}\n'.format(n_vars))
fout.write('var_names:
{}\n'.format(var_names))
for i, nm in enumerate( var_names):
fout.write(f'{nm:<.10s}
==> X[{i+1}]\n')
fout.write('\n')
X = numpy.array(X).T
RawScore = numpy.sum(X, axis=1)
meanX = numpy.mean(X, axis=0)
X = X - meanX
# Centering of X
#
# SVD of the
centered X
#
U, Sv, Vh = numpy.linalg.svd(X,
full_matrices=False)
#
# Scree plot of sigular values
#
plt.plot(numpy.arange(1, len(Sv)+1),
Sv)
plt.title('Scree Plot', fontsize= 16)
plt.ylabel('Singular value',
fontsize=14)
plt.tight_layout()
fnm_scree = s_out[:-4] + '_Scree.png'
plt.savefig(fnm_scree)
plt.show()
if numpy.sum(Vh[0]) < 0.0:
U *= -1.0
Vh *= -1.0
S = numpy.cov(X, rowvar=False)
A = Vh[0] * Sv[0] / (n_data**0.5) # Vector of Coefficients, i.e.
Lambda's
for i in range(n_vars):
print('Lambda[{0}]**2 /
Var(X[{0}]) = {1:.3f}'.format(
i+1, A[i]**2 / S[i][i]))
fout.write('Lambda[{0}]**2 /
Var(X[{0}]) = {1:.3f}\n'.format(
i+1, A[i]**2 / S[i][i]))
print()
fout.write('\n')
for i in range(n_vars):
print('{0}: Lambda[{1}] =
{2:.3f}'.format(
var_names[i], i+1, A[i]))
fout.write('{0}: Lambda[{1}]
= {2:.3f}\n'.format(
var_names[i], i+1, A[i]))
#
#
Calculation of GFI (Goodness of Fit Index)
#
sum_num = 0.0
sum_den = 0.0
for i in range(n_vars):
for j in range(n_vars):
if i
== j:
sum_den
+= S[i][i] ** 2
else:
sum_den += S[i][j] ** 2
sum_num += (S[i][j] - A[i] * A[j]) ** 2
GFI = 1.0 - sum_num / sum_den
print('\nGFI = {0:.3f}'.format(GFI))
fout.write('\nGFI =
{0:.3f}\n'.format(GFI))
#
# Calculation of omega
#
sx2 = numpy.var(X.sum(axis=1))
omega = (A.sum() ** 2) / sx2
print('\nω = {0:.3f}'.format(omega))
fout.write('\nω = {0:.3f}\n'.format(omega))
f = U[:,0] * (n_data**0.5)
f = numpy.reshape(f, (n_data,1))
A = numpy.reshape(A, (1, len(A)))
ModelScore = numpy.sum(f @ A, axis=1)
+ numpy.sum(meanX)
#
# Scatter
plot of (model prediction, raw score)
#
plt.scatter(ModelScore, RawScore)
plt.xlabel('Model Score', fontsize=14)
plt.ylabel('Raw Score', fontsize=14)
plt.plot([numpy.min(ModelScore),
numpy.max(ModelScore)],[numpy.min(ModelScore), numpy.max(ModelScore)],
label = '$Y=X$', color='k')
plt.title(f'$\\omega={omega:.3f}$',
fontsize=16)
plt.legend(fontsize=14)
plt.tight_layout()
fnm_scatter = s_out[:-4] +
'Scatter.png'
plt.savefig(fnm_scatter)
plt.show()
fout.close()
print()
print(fnm_scree + ' was saved.')
print(fnm_scatter + ' was saved.')
print(s_out + " was
saved.\n")
ベイズ分析例は、別のウェブサイトとして用意した。
McDonald, R. P. (1999). Test Theory: A unified treatment. Psychology Press.