5個の物体の重さを、1つずつ2回量る場合と、2つずつ10通り量る場合で、物体の重さが全て同じか、そうでないかF検定するとき、どちらの場合が効率的かを示しました。
コード
(1)の計画行列X1と、(2)の計画行列X2を、それぞれ0と1を反転させた計画行列X3,X4を加え、それぞれの分散と非心度λを計算して比較してみます。
# 2014 Q4(5) 2025.1.14
import numpy as np
# 設計行列 (1)~(4)
X1 = np.array([
[1, 0, 0, 0, 0],
[1, 0, 0, 0, 0],
[0, 1, 0, 0, 0],
[0, 1, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 0, 1],
[0, 0, 0, 0, 1]
])
X2 = np.array([
[1, 1, 0, 0, 0],
[1, 0, 1, 0, 0],
[1, 0, 0, 1, 0],
[1, 0, 0, 0, 1],
[0, 1, 1, 0, 0],
[0, 1, 0, 1, 0],
[0, 1, 0, 0, 1],
[0, 0, 1, 1, 0],
[0, 0, 1, 0, 1],
[0, 0, 0, 1, 1]
])
X3 = 1 - X2 # (2) の反転
X4 = 1 - X1 # (1) の反転
# 設計行列リスト
design_matrices = [X1, X2, X3, X4]
# 真の値 θ
true_theta = np.array([10, 20, 30, 40, 50])
# 結果を格納するリスト
results = []
for i, X in enumerate(design_matrices, start=1):
# 分散行列を計算
variance_matrix = np.linalg.inv(X.T @ X)
# 平均分散を計算
mean_variance = np.trace(variance_matrix) / len(variance_matrix)
# 非心度 λ を計算
theta_mean = np.mean(true_theta)
ones_vector = np.ones((X.shape[0], 1)) # X の行数に応じた列ベクトル
theta_centered = true_theta - theta_mean # θ から平均を引いたベクトル
numerator = theta_centered.T @ X.T @ X @ theta_centered
denominator = ones_vector.T @ X @ X.T @ ones_vector
non_centrality_lambda = numerator.item() # スカラーに変換
# 結果を保存
results.append((i, mean_variance, non_centrality_lambda))
# 結果の出力
for i, mean_variance, non_centrality_lambda in results:
print(f"(X{i}) の平均分散: {mean_variance:.6f}, 非心度 λ: {non_centrality_lambda:.6f}")
(X1) の平均分散: 0.500000, 非心度 λ: 2000.000000
(X2) の平均分散: 0.291667, 非心度 λ: 3000.000000
(X3) の平均分散: 0.277778, 非心度 λ: 3000.000000
(X4) の平均分散: 0.406250, 非心度 λ: 2000.000000
(2)の計画行列X2は、(1)の計画行列X1に比べ、平均分散が小さく、非心度 λが大きいため、より優れた計画行列と言えます。また、X2の0と1を反転したX3は、X2よりも分散がさらに小さく、非心度 λが変化しないことから、この中で最も優れた計画行列と考えられます。
X3 = np.array([
[0, 0, 1, 1, 1],
[0, 1, 0, 1, 1],
[0, 1, 1, 0, 1],
[0, 1, 1, 1, 0],
[1, 0, 0, 1, 1],
[1, 0, 1, 0, 1],
[1, 0, 1, 1, 0],
[1, 1, 0, 0, 1],
[1, 1, 0, 1, 0],
[1, 1, 1, 0, 0]
])