VAS値・左右瞳孔径・瞬き検出による計測（ソースコードと実行結果）

Python開発環境，ライブラリ類

ここでは、最低限の事前準備について説明する。機械学習や深層学習を行う場合は、NVIDIA CUDA、Visual Studio、Cursorなどを追加でインストールすると便利である。これらについては別ページ https://www.kkaneko.jp/cc/dev/aiassist.htmlで詳しく解説しているので、必要に応じて参照してください。

Python 3.12 のインストール

インストール済みの場合は実行不要。

管理者権限でコマンドプロンプトを起動（手順：Windowsキーまたはスタートメニュー > cmd と入力 > 右クリック > 「管理者として実行」）し、以下を実行する。管理者権限は、wingetの--scope machineオプションでシステム全体にソフトウェアをインストールするために必要である。

REM Python をシステム領域にインストール
winget install --scope machine --id Python.Python.3.12 -e --silent
REM Python のパス設定
set "PYTHON_PATH=C:\Program Files\Python312"
set "PYTHON_SCRIPTS_PATH=C:\Program Files\Python312\Scripts"
echo "%PATH%" | find /i "%PYTHON_PATH%" >nul
if errorlevel 1 setx PATH "%PATH%;%PYTHON_PATH%" /M >nul
echo "%PATH%" | find /i "%PYTHON_SCRIPTS_PATH%" >nul
if errorlevel 1 setx PATH "%PATH%;%PYTHON_SCRIPTS_PATH%" /M >nul

【関連する外部ページ】

Python の公式ページ: https://www.python.org/

AI エディタ Windsurf のインストール

Pythonプログラムの編集・実行には、AI エディタの利用を推奨する。ここでは，Windsurfのインストールを説明する。

管理者権限でコマンドプロンプトを起動（手順：Windowsキーまたはスタートメニュー > cmd と入力 > 右クリック > 「管理者として実行」）し、以下を実行して、Windsurfをシステム全体にインストールする。管理者権限は、wingetの--scope machineオプションでシステム全体にソフトウェアをインストールするために必要となる。

winget install --scope machine Codeium.Windsurf -e --silent

【関連する外部ページ】

Windsurf の公式ページ: https://windsurf.com/

必要なライブラリのインストール

コマンドプロンプトを管理者として実行（手順：Windowsキーまたはスタートメニュー > cmd と入力 > 右クリック > 「管理者として実行」）し、以下を実行する


前準備: pip install -U torch torchvision torchaudio --index-url https://download.pytorch.org/whl/cu126
pip install opencv-python scikit-learn pillow tqdm numpy

VAS値・左右瞳孔径・瞬き検出による計測プログラム


# VAS値・左右虹彩径・瞬き検出による心理生理状態計測システム
# 計測データ:
#   - VAS（Visual Analog Scale）値：0-100の主観的評価
#   - 左右虹彩径：個別に計測・正規化
#   - 虹彩径左右差：パーセンテージと状態評価
#   - 虹彩径変動係数：測定値の変動性指標
#   - 瞬き検出：タイムスタンプ付き
#   - 瞬き頻度：10秒ごとの回数と分あたり換算
# 使用技術: MediaPipe Face Mesh による虹彩輪郭計測
# 学習済みモデル: MediaPipe Face Mesh（内蔵モデル）
# 方式設計:
#   - 入力: カメラ映像、VASスライダー入力（0-100）
#   - 出力: リアルタイム表示、結果ファイル（result.txt）
#   - 処理手順: 顔画像取得→ランドマーク検出→虹彩径計測→統計分析→表示
#   - 前処理: 顔検出信頼度フィルタリング（0.5以上）
#   - 後処理: 移動平均によるノイズ除去（10フレーム）
# 注意事項: 本システムは虹彩輪郭から径を計測しており、瞳孔径の直接測定ではありません
# 前準備: pip install opencv-python mediapipe numpy matplotlib japanize-matplotlib scipy pillow

import cv2
import mediapipe as mp
import numpy as np
import tkinter as tk
from collections import deque
import time
import matplotlib.pyplot as plt
import japanize_matplotlib
from scipy import stats
from PIL import Image, ImageDraw, ImageFont
from datetime import datetime

# プログラム起動直後の解説表示
print('========== VAS値・左右虹彩径・瞬き検出による心理生理状態計測システム ==========')
print()
print('【システム概要】')
print('本システムは、Visual Analog Scale（VAS）による主観的評価と、')
print('虹彩径計測による客観的な指標を同時計測することで、')
print('意識的・無意識的な心理生理状態を包括的に評価します。')
print()
print('【Visual Analog Scale (VAS) について】')
print()
print('VASは医療・心理学分野で広く使用される主観的評価尺度です。')
print('連続的な線分上で対象者が自身の状態を示すことで、')
print('痛み、不快感、不安などの主観的体験を定量化します。')
print()
print('VASの特徴：')
print('- 0から100の連続尺度による評価')
print('- 言語化困難な感覚の数値化が可能')
print('- 時系列変化の追跡に適している')
print('- 文化や言語の影響を受けにくい')
print()
print('【VASと生理指標の関連分析】')
print()
print('本システムでは、VAS（主観的評価）と以下の客観的生理指標を同時計測し、')
print('相関関係を分析することで、以下の知見が得られます：')
print()
print('1. VAS値と虹彩径の相関')
print('   - 認知負荷や覚醒状態の変化に伴う虹彩径の変動分析')
print('   - 主観的評価と客観的指標の関係性評価')
print()
print('2. VAS値と虹彩径変動係数の関係')
print('   - 高VAS値での変動係数増大：測定環境や状態の不安定性')
print('   - 心理的動揺の計測値への影響評価')
print()
print('3. VAS値と瞬き頻度の関連')
print('   - 不快感増大時の瞬き頻度変化')
print('   - 注意資源の配分変化を反映')
print()
print('4. VAS値と左右虹彩径差の関係')
print('   - 測定条件の変化による非対称性の評価')
print('   - 計測精度の指標として活用')
print()
print('【臨床的意義】')
print()
print('- 主観的訴えの客観的指標との比較')
print('- 測定状態の安定性評価')
print('- 認知負荷の多面的評価')
print('- 測定品質の指標として活用可能')
print()
print('【計測原理と科学的根拠】')
print()
print('1. 虹彩径計測')
print('   虹彩の水平直径は約11.7±0.5mmで個人間で比較的一定です。')
print('   - 虹彩径は光の変化では変動しません（瞳孔径とは異なる）')
print('   - 個人差を正規化し、変動係数により測定安定性を定量化します')
print('   - 顔の角度や距離による見かけ上の変化を分析対象とします')
print()
print('2. 虹彩径左右差分析')
print('   健常者では左右の虹彩径はほぼ同等です（差は通常5%以内）。')
print('   10%を超える左右差は、測定条件の変化や顔の角度変化を示唆します。')
print()
print('3. 瞬き頻度分析')
print('   正常な瞬き頻度は15-20回/分です。')
print('   - 10回/分未満：過度の集中、視覚的注意の固定、疲労の蓄積')
print('   - 30回/分以上：眼精疲労、ストレス、不安状態')
print()
print('【測定値の解釈】')
print()
print('- VAS値：0（問題なし）から100（最大の不快感）の主観的評価')
print('- 虹彩径：基準値を100%とした相対値。顔の角度や距離により変動')
print('- 変動係数：虹彩径の変動性。測定環境の安定性を反映')
print('- 左右差：測定条件の対称性指標。5%以内が理想的')
print('- 瞬き頻度：覚醒度と疲労度の指標')
print()
print('【活用例】')
print('- 作業環境の客観的評価')
print('- 測定状態のモニタリング')
print('- 疲労の早期検出')
print('- 認知負荷の定量化')
print()
print('【測定上の注意事項】')
print('1. 一定の照明環境で測定してください（照明変化は計測精度に影響）')
print('2. 測定開始後10秒間は基準値設定のため安静にしてください')
print('3. コンタクトレンズや眼鏡は測定に影響しません')
print('4. 本システムは虹彩径を計測しており、瞳孔径とは異なります')
print()
print('=' * 50)
print()

# MediaPipe初期化（Python SolutionsはCPU実行が標準である）
mp_face_mesh = mp.solutions.face_mesh
face_mesh = mp_face_mesh.FaceMesh(
    static_image_mode=False,
    max_num_faces=1,
    refine_landmarks=True,
    min_detection_confidence=0.5,
    min_tracking_confidence=0.5
)
print('MediaPipe初期化完了（CPU実行）')
print()

# グローバル変数
vas_value = 50
left_iris_sizes = deque(maxlen=10)  # 虹彩径の履歴
right_iris_sizes = deque(maxlen=10)  # 虹彩径の履歴
normalized_left_sizes = deque(maxlen=10)
normalized_right_sizes = deque(maxlen=10)
baseline_left_iris_size = None  # 基準虹彩径
baseline_right_iris_size = None  # 基準虹彩径
results_data = []
last_ear = 1.0
blink_threshold = 0.2
last_blink_analysis_time = time.time()
blink_count_10sec = 0
running = True  # プログラム実行フラグ
frame_count = 0
results_log = []

# グラフ用データ
time_data = deque(maxlen=60)
vas_data_graph = deque(maxlen=60)  # 名前変更
iris_data = deque(maxlen=60)  # 虹彩径データ
cv_data = deque(maxlen=60)
lr_diff_data = deque(maxlen=60)
blink_freq_data = []
blink_freq_time = []

# 分析用データ（全期間）
all_vas_data = []
all_iris_data = []  # 虹彩径データ
all_cv_data = []
all_lr_diff_data = []
all_blink_freq_data = []
all_time_stamps = []

# 日本語フォント設定（固定）
FONT_PATH = 'C:/Windows/Fonts/meiryo.ttc'
FONT_SIZE = 20

# 虹彩ランドマークインデックスの定数化
# 左眼の虹彩輪郭ランドマーク（正しいインデックス）
LEFT_IRIS_IDX = [468, 469, 470, 471, 472]
# 右眼の虹彩輪郭ランドマーク（正しいインデックス）
RIGHT_IRIS_IDX = [473, 474, 475, 476, 477]

def create_vas_window():
    root = tk.Tk()
    root.title('VAS入力')
    root.geometry('400x100')

    def update_vas(val):
        global vas_value
        vas_value = int(float(val))

    def on_closing():
        global running
        running = False
        root.quit()

    scale = tk.Scale(root, from_=0, to=100, orient=tk.HORIZONTAL,
                     length=350, command=update_vas)
    scale.set(50)
    scale.pack(pady=20)

    label = tk.Label(root, text='0: 全く問題なし  100: 最大の不快感')
    label.pack()

    root.protocol("WM_DELETE_WINDOW", on_closing)

    return root

def compute_iris_size(indices, landmarks, image_width, image_height):
    # 楕円フィッティングによる虹彩径の推定
    points = []
    for idx in indices:
        if idx < len(landmarks):
            x = int(landmarks[idx].x * image_width)
            y = int(landmarks[idx].y * image_height)
            points.append([x, y])

    iris_size = 0
    if len(points) >= 5:
        points = np.array(points)
        # 楕円フィッティングを試みる
        try:
            ellipse = cv2.fitEllipse(points)
            # 楕円の長軸と短軸の平均を径とする
            iris_size = (ellipse[1][0] + ellipse[1][1]) / 2
        except:
            # フィッティング失敗時は点間距離の平均を使用
            distances = []
            for i in range(len(points)):
                for j in range(i+1, len(points)):
                    dist = np.linalg.norm(points[i] - points[j])
                    distances.append(dist)
            iris_size = np.mean(distances) if distances else 0
    return iris_size

def calculate_iris_size_bilateral(landmarks, image_width, image_height):
    # 左眼の虹彩輪郭ランドマーク（正しいインデックス）
    left_iris = LEFT_IRIS_IDX
    # 右眼の虹彩輪郭ランドマーク（正しいインデックス）
    right_iris = RIGHT_IRIS_IDX

    # 左眼虹彩径計算（楕円フィッティングによる改良版）
    left_iris_size = compute_iris_size(left_iris, landmarks, image_width, image_height)

    # 右眼虹彩径計算（楕円フィッティングによる改良版）
    right_iris_size = compute_iris_size(right_iris, landmarks, image_width, image_height)

    return left_iris_size, right_iris_size

def calculate_ear(landmarks):
    # 左眼のランドマーク
    left_eye = [33, 160, 158, 133, 153, 144]

    points = []
    for idx in left_eye:
        if idx < len(landmarks):
            points.append([landmarks[idx].x, landmarks[idx].y])

    if len(points) == 6:
        points = np.array(points)
        # 垂直距離
        vertical_1 = np.linalg.norm(points[1] - points[5])
        vertical_2 = np.linalg.norm(points[2] - points[4])
        # 水平距離
        horizontal = np.linalg.norm(points[0] - points[3])

        if horizontal > 0:
            ear = (vertical_1 + vertical_2) / (2.0 * horizontal)
            return ear
    return 1.0

def video_frame_processing(frame):
    global frame_count, vas_value, last_ear, blink_count_10sec
    global left_iris_sizes, right_iris_sizes, normalized_left_sizes, normalized_right_sizes
    global baseline_left_iris_size, baseline_right_iris_size
    global last_blink_analysis_time, results_data
    global time_data, vas_data_graph, iris_data, cv_data, lr_diff_data
    global blink_freq_data, blink_freq_time
    global all_vas_data, all_iris_data, all_cv_data, all_lr_diff_data
    global all_blink_freq_data, all_time_stamps

    current_time = time.time()
    frame_count += 1

    # BGR to RGB
    rgb_frame = cv2.cvtColor(frame, cv2.COLOR_BGR2RGB)
    results = face_mesh.process(rgb_frame)

    display_frame = frame.copy()
    result_text = f"フレーム{frame_count}: "

    if results.multi_face_landmarks:
        face_landmarks = results.multi_face_landmarks[0]

        # 左右虹彩径計算
        left_size, right_size = calculate_iris_size_bilateral(
            face_landmarks.landmark,
            frame.shape[1],
            frame.shape[0]
        )

        # EAR計算と瞬き検出
        current_ear = calculate_ear(face_landmarks.landmark)

        # 瞬き検出
        if last_ear > blink_threshold and current_ear <= blink_threshold:
            blink_count_10sec += 1
            result_text += f"瞬き検出 "

        last_ear = current_ear

        if left_size > 0 and right_size > 0:
            left_iris_sizes.append(left_size)
            right_iris_sizes.append(right_size)

            # ベースライン設定
            if baseline_left_iris_size is None and len(left_iris_sizes) >= 10:
                baseline_left_iris_size = np.mean(left_iris_sizes)
                baseline_right_iris_size = np.mean(right_iris_sizes)

            # 正規化と分析
            if baseline_left_iris_size and baseline_right_iris_size:
                normalized_left = (left_size / baseline_left_iris_size) * 100
                normalized_right = (right_size / baseline_right_iris_size) * 100
                normalized_left_sizes.append(normalized_left)
                normalized_right_sizes.append(normalized_right)

                # 平均虹彩径と変動係数
                avg_normalized = (normalized_left + normalized_right) / 2
                if len(normalized_left_sizes) >= 5:
                    all_normalized = list(normalized_left_sizes) + list(normalized_right_sizes)
                    cv_value = np.std(all_normalized) / np.mean(all_normalized) * 100
                else:
                    cv_value = 0

                # 左右差計算（絶対値）
                if (normalized_left + normalized_right) > 0:
                    lr_diff = abs(normalized_left - normalized_right) / ((normalized_left + normalized_right) / 2) * 100
                    if lr_diff <= 5:
                        lr_status = '理想的'
                    elif lr_diff <= 10:
                        lr_status = '軽度非対称'
                    else:
                        lr_status = '要確認'
                else:
                    lr_diff = 0
                    lr_status = '計測中'

                # ランドマークをプロット
                # 左眼虹彩ランドマーク（青色）
                left_iris_landmarks = LEFT_IRIS_IDX
                for idx in left_iris_landmarks:
                    if idx < len(face_landmarks.landmark):
                        x = int(face_landmarks.landmark[idx].x * frame.shape[1])
                        y = int(face_landmarks.landmark[idx].y * frame.shape[0])
                        cv2.circle(display_frame, (x, y), 8, (255, 0, 0), -1)

                # 右眼虹彩ランドマーク（赤色）
                right_iris_landmarks = RIGHT_IRIS_IDX
                for idx in right_iris_landmarks:
                    if idx < len(face_landmarks.landmark):
                        x = int(face_landmarks.landmark[idx].x * frame.shape[1])
                        y = int(face_landmarks.landmark[idx].y * frame.shape[0])
                        cv2.circle(display_frame, (x, y), 8, (0, 0, 255), -1)

                # 瞬き検出用ランドマーク（緑色）
                blink_landmarks = [33, 160, 158, 133, 153, 144]
                for idx in blink_landmarks:
                    if idx < len(face_landmarks.landmark):
                        x = int(face_landmarks.landmark[idx].x * frame.shape[1])
                        y = int(face_landmarks.landmark[idx].y * frame.shape[0])
                        cv2.circle(display_frame, (x, y), 5, (0, 255, 0), -1)

                # 画面表示（Pillow+OpenCV）
                img_pil = Image.fromarray(cv2.cvtColor(display_frame, cv2.COLOR_BGR2RGB))
                draw = ImageDraw.Draw(img_pil)
                font = ImageFont.truetype(FONT_PATH, FONT_SIZE)
                draw.text((10, 30), f'VAS: {vas_value}', font=font, fill=(0, 255, 0))
                draw.text((10, 60), f'虹彩径: {avg_normalized:.1f}%', font=font, fill=(0, 255, 0))
                draw.text((10, 90), f'変動係数: {cv_value:.1f}%', font=font, fill=(0, 255, 0))
                draw.text((10, 120), f'左右差: {lr_diff:.1f}% ({lr_status})', font=font, fill=(0, 255, 0))
                draw.text((10, 150), f'青点:左眼虹彩 赤点:右眼虹彩 緑点:瞬き検出', font=font, fill=(255, 255, 255))
                display_frame = cv2.cvtColor(np.array(img_pil), cv2.COLOR_RGB2BGR)

                result_text += f"VAS:{vas_value} 虹彩径:{avg_normalized:.1f}% CV:{cv_value:.1f}% 左右差:{lr_diff:.1f}%"

                # データ保存
                time_data.append(current_time)
                vas_data_graph.append(vas_value)
                iris_data.append(avg_normalized)
                cv_data.append(cv_value)
                lr_diff_data.append(lr_diff)

                all_time_stamps.append(current_time)
                all_vas_data.append(vas_value)
                all_iris_data.append(avg_normalized)
                all_cv_data.append(cv_value)
                all_lr_diff_data.append(lr_diff)

                # 10秒ごとの瞬き頻度分析
                elapsed_time = current_time - last_blink_analysis_time
                if elapsed_time >= 10.0:
                    # 実際の経過時間に基づいて分あたりに換算
                    blink_per_minute = (blink_count_10sec / elapsed_time) * 60

                    blink_freq_data.append(blink_per_minute)
                    blink_freq_time.append(current_time)
                    all_blink_freq_data.append(blink_per_minute)

                    result_text += f" 瞬き頻度:{blink_per_minute:.0f}回/分"

                    blink_count_10sec = 0
                    last_blink_analysis_time = current_time
    else:
        result_text += "顔検出失敗"

    # グラフ更新
    update_graphs()

    return display_frame, result_text, current_time

# グラフ初期化
plt.ion()
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(12, 8))
fig.suptitle('心理生理状態リアルタイムモニタリング', fontsize=16)

# グラフ1: 主観評価と虹彩径
ax1.set_xlabel('時間（秒）')
ax1.set_ylabel('値')
ax1.set_title('主観評価と虹彩径の時系列変化')
ax1.set_ylim(0, 150)
line1_vas, = ax1.plot([], [], 'b-', label='VAS値')
line1_iris, = ax1.plot([], [], 'r-', label='虹彩径（%）')
ax1.legend()
ax1.text(0.5, -0.15, '青線：主観的不快感、赤線：虹彩径の相対値',
         transform=ax1.transAxes, ha='center', fontsize=10)

# グラフ2: 変動指標
ax2.set_xlabel('時間（秒）')
ax2.set_ylabel('値')
ax2.set_title('変動指標')
ax2.set_ylim(0, 50)
line2_cv, = ax2.plot([], [], 'g-', label='変動係数（%）')
ax2_twin = ax2.twinx()
ax2_twin.set_ylabel('瞬き頻度（回/分）')
ax2_twin.set_ylim(0, 60)
line2_blink, = ax2_twin.plot([], [], 'mo', markersize=8, label='瞬き頻度')
ax2.legend(loc='upper left')
ax2_twin.legend(loc='upper right')
ax2.text(0.5, -0.15, '緑線：虹彩径の変動性、紫点：10秒ごとの瞬き頻度',
         transform=ax2.transAxes, ha='center', fontsize=10)

# グラフ3: 左右対称性
ax3.set_xlabel('時間（秒）')
ax3.set_ylabel('左右差（%）')
ax3.set_title('虹彩径左右差の時系列変化')
ax3.set_ylim(0, 30)
line3, = ax3.plot([], [], 'orange')
ax3.axhline(y=5, color='green', linestyle='--', alpha=0.5, label='理想範囲')
ax3.axhline(y=10, color='red', linestyle='--', alpha=0.5, label='要確認')
ax3.legend()
ax3.text(0.5, -0.15, '5%以内が理想的、10%超は測定条件要確認',
         transform=ax3.transAxes, ha='center', fontsize=10)

# グラフ4: 瞬き頻度トレンド
ax4.set_xlabel('時間（秒）')
ax4.set_ylabel('瞬き頻度（回/分）')
ax4.set_title('瞬き頻度の推移')
ax4.set_ylim(0, 60)
line4, = ax4.plot([], [], 'bo-', markersize=6)
ax4.axhspan(15, 20, alpha=0.2, color='green', label='正常範囲')
ax4.legend()
ax4.text(0.5, -0.15, '15-20回/分が正常、10回未満は疲労、30回超はストレス',
         transform=ax4.transAxes, ha='center', fontsize=10)

plt.tight_layout()

def update_graphs():
    try:
        if len(time_data) > 0:
            time_array = np.array(time_data) - time_data[0]

            # グラフ1更新
            line1_vas.set_data(time_array, vas_data_graph)
            line1_iris.set_data(time_array, iris_data)
            ax1.set_xlim(max(0, time_array[-1] - 60), time_array[-1] + 1)

            # グラフ2更新
            line2_cv.set_data(time_array, cv_data)
            if len(blink_freq_time) > 0:
                blink_time_array = np.array(blink_freq_time) - time_data[0]
                line2_blink.set_data(blink_time_array, blink_freq_data)
            ax2.set_xlim(max(0, time_array[-1] - 60), time_array[-1] + 1)
            ax2_twin.set_xlim(max(0, time_array[-1] - 60), time_array[-1] + 1)

            # グラフ3更新
            line3.set_data(time_array, lr_diff_data)
            ax3.set_xlim(max(0, time_array[-1] - 60), time_array[-1] + 1)

            # グラフ4更新
            if len(blink_freq_time) > 0:
                blink_time_array = np.array(blink_freq_time) - time_data[0]
                line4.set_data(blink_time_array, blink_freq_data)
                ax4.set_xlim(max(0, time_array[-1] - 60), time_array[-1] + 1)

        plt.draw()
        plt.pause(0.001)
    except:
        pass

def calculate_correlations():
    print('\n===== データ分析結果 =====')
    print(f'測定時間: {len(all_time_stamps)}秒')
    print(f'データ点数: {len(all_vas_data)}個')
    print()
    print('【相関の解釈基準】')
    print('|r| ≥ 0.7: 強い相関')
    print('0.4 ≤ |r| < 0.7: 中程度の相関')
    print('0.2 ≤ |r| < 0.4: 弱い相関')
    print('|r| < 0.2: ほぼ相関なし')
    print('p < 0.05: 統計的に有意')
    print()
    print('【相関分析結果】')

    correlation_results = []
    correlation_results.append('===== データ分析結果 =====')
    correlation_results.append(f'測定時間: {len(all_time_stamps)}秒')
    correlation_results.append(f'データ点数: {len(all_vas_data)}個')
    correlation_results.append('')
    correlation_results.append('【相関の解釈基準】')
    correlation_results.append('|r| ≥ 0.7: 強い相関')
    correlation_results.append('0.4 ≤ |r| < 0.7: 中程度の相関')
    correlation_results.append('0.2 ≤ |r| < 0.4: 弱い相関')
    correlation_results.append('|r| < 0.2: ほぼ相関なし')
    correlation_results.append('p < 0.05: 統計的に有意')
    correlation_results.append('')
    correlation_results.append('【相関分析結果】')

    # 相関分析の実行
    if len(all_vas_data) > 2 and len(all_iris_data) > 2:
        r1, p1 = stats.pearsonr(all_vas_data, all_iris_data)
        interpretation1 = interpret_correlation(r1)
        print(f'1. VAS値 × 平均虹彩径')
        print(f'   相関係数: {r1:.3f} (p値: {p1:.3f})')
        print(f'   解釈: {interpretation1}')
        print()

        correlation_results.append(f'1. VAS値 × 平均虹彩径')
        correlation_results.append(f'   相関係数: {r1:.3f} (p値: {p1:.3f})')
        correlation_results.append(f'   解釈: {interpretation1}')
        correlation_results.append('')

    if len(all_vas_data) > 2 and len(all_cv_data) > 2:
        r2, p2 = stats.pearsonr(all_vas_data, all_cv_data)
        interpretation2 = interpret_correlation(r2)
        print(f'2. VAS値 × 虹彩径変動係数')
        print(f'   相関係数: {r2:.3f} (p値: {p2:.3f})')
        print(f'   解釈: {interpretation2}')
        print()

        correlation_results.append(f'2. VAS値 × 虹彩径変動係数')
        correlation_results.append(f'   相関係数: {r2:.3f} (p値: {p2:.3f})')
        correlation_results.append(f'   解釈: {interpretation2}')
        correlation_results.append('')

    if len(blink_freq_time) > 2 and len(all_vas_data) > 2:
        vas_at_blink_times = np.interp(blink_freq_time, all_time_stamps, all_vas_data)
        min_len = min(len(vas_at_blink_times), len(all_blink_freq_data))
        if min_len > 2:
            r3, p3 = stats.pearsonr(vas_at_blink_times[:min_len], all_blink_freq_data[:min_len])
            interpretation3 = interpret_correlation(r3)
            print(f'3. VAS値 × 瞬き頻度')
            print(f'   相関係数: {r3:.3f} (p値: {p3:.3f})')
            print(f'   解釈: {interpretation3}')
            print()

            correlation_results.append(f'3. VAS値 × 瞬き頻度')
            correlation_results.append(f'   相関係数: {r3:.3f} (p値: {p3:.3f})')
            correlation_results.append(f'   解釈: {interpretation3}')
            correlation_results.append('')

    if len(blink_freq_time) > 2 and len(all_cv_data) > 2:
        cv_at_blink_times = np.interp(blink_freq_time, all_time_stamps, all_cv_data)
        min_len = min(len(cv_at_blink_times), len(all_blink_freq_data))
        if min_len > 2:
            r4, p4 = stats.pearsonr(cv_at_blink_times[:min_len], all_blink_freq_data[:min_len])
            interpretation4 = interpret_correlation(r4)
            print(f'4. 虹彩径変動係数 × 瞬き頻度')
            print(f'   相関係数: {r4:.3f} (p値: {p4:.3f})')
            print(f'   解釈: {interpretation4}')
            print()

            correlation_results.append(f'4. 虹彩径変動係数 × 瞬き頻度')
            correlation_results.append(f'   相関係数: {r4:.3f} (p値: {p4:.3f})')
            correlation_results.append(f'   解釈: {interpretation4}')
            correlation_results.append('')

    if len(all_vas_data) > 2 and len(all_lr_diff_data) > 2:
        r5, p5 = stats.pearsonr(all_vas_data, all_lr_diff_data)
        interpretation5 = interpret_correlation(r5)
        print(f'5. VAS値 × 虹彩径左右差')
        print(f'   相関係数: {r5:.3f} (p値: {p5:.3f})')
        print(f'   解釈: {interpretation5}')
        print()

        correlation_results.append(f'5. VAS値 × 虹彩径左右差')
        correlation_results.append(f'   相関係数: {r5:.3f} (p値: {p5:.3f})')
        correlation_results.append(f'   解釈: {interpretation5}')
        correlation_results.append('')

    # 統計サマリー
    print('【統計サマリー】')
    correlation_results.append('【統計サマリー】')

    if len(all_vas_data) > 0:
        print(f'- VAS値: 平均 {np.mean(all_vas_data):.1f}, 標準偏差 {np.std(all_vas_data):.1f}, 範囲 {min(all_vas_data)}-{max(all_vas_data)}')
        correlation_results.append(f'- VAS値: 平均 {np.mean(all_vas_data):.1f}, 標準偏差 {np.std(all_vas_data):.1f}, 範囲 {min(all_vas_data)}-{max(all_vas_data)}')

    if len(all_iris_data) > 0:
        print(f'- 虹彩径: 平均 {np.mean(all_iris_data):.1f}%, 標準偏差 {np.std(all_iris_data):.1f}%, 範囲 {min(all_iris_data):.1f}-{max(all_iris_data):.1f}%')
        correlation_results.append(f'- 虹彩径: 平均 {np.mean(all_iris_data):.1f}%, 標準偏差 {np.std(all_iris_data):.1f}%, 範囲 {min(all_iris_data):.1f}-{max(all_iris_data):.1f}%')

    if len(all_cv_data) > 0:
        print(f'- 変動係数: 平均 {np.mean(all_cv_data):.1f}%, 標準偏差 {np.std(all_cv_data):.1f}%, 範囲 {min(all_cv_data):.1f}-{max(all_cv_data):.1f}%')
        correlation_results.append(f'- 変動係数: 平均 {np.mean(all_cv_data):.1f}%, 標準偏差 {np.std(all_cv_data):.1f}%, 範囲 {min(all_cv_data):.1f}-{max(all_cv_data):.1f}%')

    if len(all_blink_freq_data) > 0:
        print(f'- 瞬き頻度: 平均 {np.mean(all_blink_freq_data):.1f}回/分, 標準偏差 {np.std(all_blink_freq_data):.1f}, 範囲 {min(all_blink_freq_data):.0f}-{max(all_blink_freq_data):.0f}回/分')
        correlation_results.append(f'- 瞬き頻度: 平均 {np.mean(all_blink_freq_data):.1f}回/分, 標準偏差 {np.std(all_blink_freq_data):.1f}, 範囲 {min(all_blink_freq_data):.0f}-{max(all_blink_freq_data):.0f}回/分')

    if len(all_lr_diff_data) > 0:
        print(f'- 左右差: 平均 {np.mean(all_lr_diff_data):.1f}%, 標準偏差 {np.std(all_lr_diff_data):.1f}%, 最大 {max(all_lr_diff_data):.1f}%')
        correlation_results.append(f'- 左右差: 平均 {np.mean(all_lr_diff_data):.1f}%, 標準偏差 {np.std(all_lr_diff_data):.1f}%, 最大 {max(all_lr_diff_data):.1f}%')

    return correlation_results

def interpret_correlation(r):
    abs_r = abs(r)
    if abs_r >= 0.7:
        strength = '強い'
    elif abs_r >= 0.4:
        strength = '中程度の'
    elif abs_r >= 0.2:
        strength = '弱い'
    else:
        return 'ほぼ相関なし'

    direction = '正の' if r > 0 else '負の'
    return f'{direction}{strength}相関あり'

def print_usage():
    print('操作方法:')
    print('- VASスライダーで主観的な状態を入力してください（0-100）')
    print('- qキーで終了します')
    print('- 画面表示: 青点=左眼虹彩、赤点=右眼虹彩、緑点=瞬き検出用ランドマーク')
    print()

print_usage()

# VASウィンドウ作成
vas_window = create_vas_window()

# カメラ初期化
cap = cv2.VideoCapture(0, cv2.CAP_DSHOW)
if not cap.isOpened():
    cap = cv2.VideoCapture(0)
cap.set(cv2.CAP_PROP_BUFFERSIZE, 1)

if not cap.isOpened():
    print('カメラの初期化に失敗しました')
    exit()

# 計測開始時の説明表示
print('===== 計測システムの詳細 =====')
print()
print('【計測データ】')
print('本システムは以下の6種類のデータを同時計測します：')
print()
print('1. VAS（Visual Analog Scale）値')
print('   - 0-100の連続尺度による主観的評価')
print('   - スライダーで随時調整可能')
print()
print('2. 左右虹彩径')
print('   - MediaPipe Face Meshの虹彩ランドマーク（左眼:468-472、右眼:473-477）を使用')
print('   - 虹彩輪郭から径を計測')
print('   - 個人差補正のため初期10フレームで正規化')
print()
print('3. 虹彩径左右差')
print('   - (|左眼径-右眼径|)/平均径×100で算出')
print('   - 5%以内:理想的、5-10%:軽度非対称、10%超:要確認')
print()
print('4. 虹彩径変動係数（CV）')
print('   - 標準偏差/平均×100で算出')
print('   - 測定値の変動性を定量化')
print()
print('5. 瞬き検出')
print('   - Eye Aspect Ratio（EAR）による検出')
print('   - 発生時刻をミリ秒単位で記録')
print()
print('6. 瞬き頻度')
print('   - 10秒間隔で集計し分あたりに換算')
print('   - 15-20回/分が正常範囲')
print()
print('【計測原理】')
print('虹彩径の測定値の変動を分析し、')
print('VASとの相関分析により、主観と客観の関係を評価します。')
print()
print('注意：本システムは虹彩輪郭から径を計測しており、')
print('瞳孔径の直接測定ではありません。')
print()
print('計測を開始します...')
print('=====================================')
print()

print('\n=== 動画処理開始 ===')
print_usage()

try:
    while running:
        try:
            vas_window.update_idletasks()
            vas_window.update()
        except:
            running = False
            break

        ret, frame = cap.read()
        if not ret:
            break

        MAIN_FUNC_DESC = "VAS・虹彩径・瞬き計測"
        processed_frame, result, current_time = video_frame_processing(frame)
        cv2.imshow(MAIN_FUNC_DESC, processed_frame)

        print(datetime.fromtimestamp(current_time).strftime("%Y-%m-%d %H:%M:%S.%f")[:-3], result)
        results_log.append(result)

        if cv2.waitKey(1) & 0xFF == ord('q'):
            running = False
            break

finally:
    print('\n=== プログラム終了 ===')
    cap.release()
    cv2.destroyAllWindows()
    try:
        vas_window.destroy()
    except:
        pass
    plt.close('all')

    # 相関分析
    correlation_results = calculate_correlations()

    # 結果保存
    if results_log:
        with open('result.txt', 'w', encoding='utf-8') as f:
            f.write('=== 結果 ===\n')
            f.write(f'処理フレーム数: {frame_count}\n')
            f.write(f'使用デバイス: CPU\n')
            f.write('\n')
            f.write('\n'.join(results_log))
            f.write('\n\n')
            for line in correlation_results:
                f.write(line + '\n')
        print(f'\n処理結果をresult.txtに保存しました')