Files

593 lines
21 KiB
Python
Raw Permalink Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
#!/usr/bin/env python3
# =============================================================================
# FeatureService 3D — Многоканальная радиомика Brain Tumour
# РЕЗУЛЬТАТЫ СОХРАНЯЮТСЯ РЯДОМ С ДАННЫМИ
# =============================================================================
import os
import sys
import time
import logging
import json
import numpy as np
import pandas as pd
from pathlib import Path
from tqdm import tqdm
import matplotlib.pyplot as plt
import seaborn as sns
# PyRadiomics
import radiomics
from radiomics import featureextractor
import SimpleITK as sitk
import nibabel as nib
# =============================================================================
# КОНФИГУРАЦИЯ
# =============================================================================
# --- Путь к данным ---
DATA_ROOT = Path(r"D:\Radiomics_Data\Experiments\experiment-17-06-2026\Task01_BrainTumour\Task01_BrainTumour")
# --- Пути к данным внутри DATA_ROOT ---
PATH_IMAGES_TRAIN = DATA_ROOT / "imagesTr"
PATH_IMAGES_TEST = DATA_ROOT / "imagesTs"
PATH_LABELS_TRAIN = DATA_ROOT / "labelsTr"
# --- !!! НОВЫЕ ПУТИ ДЛЯ РЕЗУЛЬТАТОВ (ВНУТРИ DATA_ROOT) !!! ---
OUTPUT_ROOT = DATA_ROOT / "radiomics_output" # ← Теперь рядом с данными!
PATH_FEATURES = OUTPUT_ROOT / "features"
PATH_VISUALIZATION = OUTPUT_ROOT / "visualization"
PATH_LOGS = OUTPUT_ROOT / "logs"
PATH_MODELS = OUTPUT_ROOT / "models"
# --- Проект (только для справки) ---
PROJECT_ROOT = Path(__file__).parent.resolve()
FILE_EXTENSION = ".nii.gz"
MAX_PATIENTS = 2
# Названия модальностей (каналов)
CHANNEL_NAMES = ['T1', 'T1c', 'T2', 'FLAIR']
# =============================================================================
# СОЗДАНИЕ ПАПОК
# =============================================================================
def ensure_directories():
"""Создаёт все необходимые директории внутри папки с данными"""
print("\n📁 СОЗДАНИЕ ПАПОК ДЛЯ РЕЗУЛЬТАТОВ:")
for path in [PATH_FEATURES, PATH_VISUALIZATION, PATH_LOGS, PATH_MODELS]:
path.mkdir(parents=True, exist_ok=True)
print(f" ✅ {path}")
print()
ensure_directories()
# =============================================================================
# НАСТРОЙКА ЛОГГИРОВАНИЯ
# =============================================================================
logging.basicConfig(
level=logging.INFO,
format='%(asctime)s - %(levelname)s - %(message)s',
handlers=[
logging.FileHandler(PATH_LOGS / 'radiomics_3d.log', encoding='utf-8'),
logging.StreamHandler()
]
)
logger = logging.getLogger(__name__)
# =============================================================================
# НАСТРОЙКИ RADIOMICS
# =============================================================================
RADIOMICS_SETTINGS = {
'minimumROIDimensions': 2,
'minimumROISize': None,
'normalize': False,
'normalizeScale': 1,
'removeOutliers': None,
'resampledPixelSpacing': None,
'interpolator': 'sitkBSpline',
'preCrop': True,
'padDistance': 5,
'distances': [1],
'force2D': False,
'force2Ddimension': 0,
'resegmentRange': None,
'label': 1
}
FEATURE_CLASSES = [
'firstorder',
'shape',
'glcm',
'glrlm',
'glszm',
'gldm',
'ngtdm'
]
ENABLED_IMAGE_TYPES = {
'Original': {},
'Wavelet': {},
'Square': {},
'SquareRoot': {},
'Logarithm': {},
'Exponential': {},
'Gradient': {}
}
# =============================================================================
# ВСПОМОГАТЕЛЬНЫЕ ФУНКЦИИ
# =============================================================================
def print_structure():
"""Выводит структуру проекта и данных"""
print("\n" + "=" * 70)
print("🧠 3D РАДИОМИКА — BRAIN TUMOUR (МНОГОКАНАЛЬНАЯ)")
print("=" * 70)
print(f" 📂 Папка с данными: {DATA_ROOT}")
print(f" 📂 Изображения: {PATH_IMAGES_TRAIN}")
print(f" 📂 Маски: {PATH_LABELS_TRAIN}")
print(f" 📂 РЕЗУЛЬТАТЫ: {OUTPUT_ROOT}")
print(f" ├── features/ {PATH_FEATURES}")
print(f" ├── visualization/{PATH_VISUALIZATION}")
print(f" └── logs/ {PATH_LOGS}")
print("=" * 70 + "\n")
def check_data():
"""Проверяет наличие данных"""
print("🔍 ПРОВЕРКА ДАННЫХ:")
print("-" * 50)
for name, path in [
("Обучающие изображения", PATH_IMAGES_TRAIN),
("Тестовые изображения", PATH_IMAGES_TEST),
("Маски", PATH_LABELS_TRAIN)
]:
if path.exists():
files = [f for f in path.glob(f'*{FILE_EXTENSION}')
if not f.name.startswith('._')]
print(f" ✅ {name}: {len(files)} файлов")
else:
print(f" ❌ {name}: папка не найдена")
# Проверяем, создалась ли папка для результатов
if OUTPUT_ROOT.exists():
print(f" ✅ Папка результатов: {OUTPUT_ROOT}")
else:
print(f" ❌ Папка результатов не создана")
print("-" * 50 + "\n")
def is_valid_nifti_file(filepath):
name = filepath.name
if name.startswith('._'):
return False
if not name.endswith('.nii.gz'):
return False
if filepath.stat().st_size < 1000:
return False
return True
def load_nifti_channels(filepath):
"""Загружает NIfTI и возвращает все каналы"""
try:
if not Path(filepath).exists():
return None, None
img = nib.load(str(filepath))
data = img.get_fdata()
if data.ndim == 4:
num_channels = data.shape[3]
logger.info(f" 📊 Найдено каналов: {num_channels}")
channels = []
for i in range(num_channels):
channel_data = data[:, :, :, i]
if np.isnan(channel_data).any():
channel_data = np.nan_to_num(channel_data)
channels.append(channel_data)
return channels, num_channels
elif data.ndim == 3:
logger.info(f" 📊 Одноканальное изображение")
if np.isnan(data).any():
data = np.nan_to_num(data)
return [data], 1
else:
logger.error(f"❌ Неподдерживаемая размерность: {data.ndim}")
return None, None
except Exception as e:
logger.error(f"❌ Ошибка загрузки {filepath}: {e}")
return None, None
def create_sitk_image_from_array(data):
"""Создаёт SimpleITK изображение из массива"""
if np.isnan(data).any():
data = np.nan_to_num(data)
return sitk.GetImageFromArray(data)
# =============================================================================
# ОСНОВНОЙ КЛАСС
# =============================================================================
class MultiChannelRadiomicsExtractor:
def __init__(self, settings=None, feature_classes=None, image_types=None):
self.settings = settings or RADIOMICS_SETTINGS
self.feature_classes = feature_classes or FEATURE_CLASSES
self.image_types = image_types or ENABLED_IMAGE_TYPES
self.channel_names = CHANNEL_NAMES
logger.info("🔧 Инициализация многоканального экстрактора...")
try:
self.extractor = featureextractor.RadiomicsFeatureExtractor()
for key, value in self.settings.items():
self.extractor.settings[key] = value
for feature_class in self.feature_classes:
try:
self.extractor.enableFeatureClassByName(feature_class)
logger.info(f" ✅ Включён класс: {feature_class}")
except Exception as e:
logger.warning(f" ⚠️ Не удалось включить {feature_class}: {e}")
self.extractor.enableImageTypes(**self.image_types)
logger.info("✅ Экстрактор инициализирован")
except Exception as e:
logger.error(f"❌ Ошибка инициализации: {e}")
raise
def extract_features_for_channel(self, image_data, mask_data, patient_id, channel_name):
"""Извлекает признаки для одного канала"""
try:
image_sitk = create_sitk_image_from_array(image_data)
mask_sitk = create_sitk_image_from_array(mask_data)
mask_array = sitk.GetArrayFromImage(mask_sitk)
if mask_array.sum() == 0:
logger.warning(f" ⚠️ Пустая маска для {channel_name}")
return None
result = self.extractor.execute(image_sitk, mask_sitk)
features = {}
for key, value in result.items():
if isinstance(value, (int, float, str, bool, np.number)):
new_key = f"{channel_name}_{key}"
features[new_key] = value
elif isinstance(value, np.ndarray):
new_key = f"{channel_name}_{key}"
features[new_key] = value.tolist()
return features
except Exception as e:
logger.error(f" ❌ Ошибка для {channel_name}: {e}")
return None
def extract_all_channels(self, image_path, mask_path, patient_id):
"""Извлекает признаки для всех каналов"""
logger.info(f"🔄 Обработка: {patient_id}")
channels, num_channels = load_nifti_channels(image_path)
if channels is None:
logger.error(f"❌ Не удалось загрузить изображение")
return None
mask_channels, _ = load_nifti_channels(mask_path)
if mask_channels is None:
logger.error(f"❌ Не удалось загрузить маску")
return None
mask_data = mask_channels[0]
channel_names = self.channel_names[:len(channels)]
if len(channel_names) < len(channels):
channel_names = [f"Channel_{i}" for i in range(len(channels))]
all_features = {
'patient_id': patient_id,
'image_path': str(image_path),
'mask_path': str(mask_path),
'num_channels': len(channels)
}
start_time = time.time()
for i, (channel_data, ch_name) in enumerate(zip(channels, channel_names)):
logger.info(f" 📊 Канал {i + 1}: {ch_name}")
features = self.extract_features_for_channel(
channel_data, mask_data, patient_id, ch_name
)
if features:
all_features.update(features)
else:
logger.warning(f" ⚠️ Не удалось извлечь признаки для {ch_name}")
elapsed = time.time() - start_time
logger.info(f"✅ {patient_id}: {len(all_features)} признаков ({elapsed:.2f}с)")
return all_features
def process_patient(self, image_path, mask_path, output_dir=None):
if output_dir is None:
output_dir = PATH_FEATURES
patient_id = Path(image_path).stem.replace('.nii', '')
features = self.extract_all_channels(image_path, mask_path, patient_id)
if features is None:
return None
output_path = Path(output_dir) / f"{patient_id}_multichannel_features.csv"
try:
output_path.parent.mkdir(parents=True, exist_ok=True)
pd.DataFrame([features]).to_csv(output_path, index=False, encoding='utf-8')
logger.info(f"💾 Сохранено: {output_path}")
return str(output_path)
except Exception as e:
logger.error(f"❌ Ошибка сохранения: {e}")
return None
def process_all_patients(self, images_dir, masks_dir, output_dir=None):
if output_dir is None:
output_dir = PATH_FEATURES
if not Path(images_dir).exists():
logger.error(f"❌ Папка не найдена: {images_dir}")
return []
if not Path(masks_dir).exists():
logger.error(f"❌ Папка не найдена: {masks_dir}")
return []
all_files = list(Path(images_dir).glob(f'*{FILE_EXTENSION}'))
image_files = [f for f in all_files if is_valid_nifti_file(f)]
if MAX_PATIENTS:
image_files = image_files[:MAX_PATIENTS]
logger.info(f"📁 Найдено пациентов: {len(image_files)}")
if not image_files:
logger.warning("⚠️ Нет файлов")
return []
results = []
for img_path in tqdm(image_files, desc="Обработка пациентов"):
patient_id = img_path.stem.replace('.nii', '')
mask_path = Path(masks_dir) / f"{patient_id}{FILE_EXTENSION}"
if not mask_path.exists():
logger.warning(f"⚠️ Нет маски для {patient_id}")
continue
if not is_valid_nifti_file(mask_path):
logger.warning(f"⚠️ Маска системный файл")
continue
result = self.process_patient(img_path, mask_path, output_dir)
if result:
results.append(result)
logger.info(f"✅ Обработано: {len(results)} пациентов")
return results
# =============================================================================
# ВИЗУАЛИЗАЦИЯ
# =============================================================================
def visualize_features(features_df, output_dir=None):
if output_dir is None:
output_dir = PATH_VISUALIZATION
output_dir = Path(output_dir)
output_dir.mkdir(parents=True, exist_ok=True)
if features_df is None or len(features_df) == 0:
logger.warning("⚠️ Нет данных для визуализации")
return
key_features = [
'original_shape_Volume',
'original_shape_SurfaceArea',
'original_shape_Sphericity',
'original_firstorder_Mean',
'original_firstorder_Entropy',
'original_firstorder_Skewness',
'original_glcm_Contrast',
'original_glcm_Correlation',
'original_glrlm_RunEntropy',
'original_glszm_ZoneEntropy'
]
available_features = [col for col in features_df.columns
if any(kf in col for kf in key_features)]
if not available_features:
logger.warning("⚠️ Нет доступных признаков для визуализации")
return
print("\n" + "=" * 60)
print("📊 ВИЗУАЛИЗАЦИЯ РЕЗУЛЬТАТОВ")
print("=" * 60)
print(f" 📂 Сохранение в: {output_dir}")
# Распределения
fig, axes = plt.subplots(3, 4, figsize=(16, 12))
axes = axes.flatten()
for i, feature in enumerate(available_features[:12]):
if i < len(axes):
sns.histplot(features_df[feature].dropna(), kde=True, ax=axes[i])
axes[i].set_title(feature[:40], fontsize=9)
axes[i].set_xlabel('')
for j in range(i + 1, len(axes)):
axes[j].set_visible(False)
plt.tight_layout()
plt.savefig(output_dir / 'feature_distributions.png', dpi=300, bbox_inches='tight')
plt.close()
print(f" ✅ Распределения: {output_dir / 'feature_distributions.png'}")
# Корреляции
top_features = available_features[:10]
if len(top_features) > 1:
plt.figure(figsize=(12, 10))
corr_matrix = features_df[top_features].corr()
sns.heatmap(corr_matrix, annot=True, fmt='.2f', cmap='coolwarm',
square=True, annot_kws={'size': 8})
plt.title('Корреляционная матрица признаков')
plt.tight_layout()
plt.savefig(output_dir / 'correlation_matrix.png', dpi=300, bbox_inches='tight')
plt.close()
print(f" ✅ Корреляции: {output_dir / 'correlation_matrix.png'}")
# Сравнение каналов
channel_cols = [col for col in features_df.columns
if any(f'_{ch}_' in col for ch in ['T1', 'T1c', 'T2', 'FLAIR'])]
if channel_cols:
feature_types = set()
for col in channel_cols:
for ft in ['firstorder_Mean', 'glcm_Contrast', 'glcm_Entropy']:
if ft in col:
feature_types.add(ft)
for ft in feature_types:
ft_cols = [col for col in channel_cols if ft in col]
if len(ft_cols) > 1:
plt.figure(figsize=(10, 6))
data_to_plot = [features_df[col].dropna() for col in ft_cols]
channels = [col.split('_')[0] for col in ft_cols]
plt.boxplot(data_to_plot, labels=channels)
plt.title(f'Сравнение каналов: {ft}')
plt.ylabel(ft)
plt.tight_layout()
plt.savefig(output_dir / f'channel_comparison_{ft}.png', dpi=300, bbox_inches='tight')
plt.close()
print(f" ✅ Сравнение каналов: {output_dir / f'channel_comparison_{ft}.png'}")
# =============================================================================
# ОБЪЕДИНЕНИЕ РЕЗУЛЬТАТОВ
# =============================================================================
def merge_features(features_dir=None, output_path=None):
if features_dir is None:
features_dir = PATH_FEATURES
if output_path is None:
output_path = features_dir / "combined_multichannel_features.csv"
features_dir = Path(features_dir)
csv_files = list(features_dir.glob("*_multichannel_features.csv"))
if not csv_files:
logger.warning("⚠️ Нет CSV файлов для объединения")
return None
logger.info(f"🔗 Объединяю {len(csv_files)} файлов...")
all_dfs = []
for csv_file in tqdm(csv_files, desc="Объединение"):
try:
df = pd.read_csv(csv_file)
all_dfs.append(df)
except Exception as e:
logger.error(f"❌ Ошибка чтения {csv_file}: {e}")
if not all_dfs:
return None
combined_df = pd.concat(all_dfs, ignore_index=True)
combined_df.to_csv(output_path, index=False, encoding='utf-8')
logger.info(f"💾 Объединённый файл: {output_path}")
logger.info(f" Пациентов: {len(combined_df)}")
logger.info(f" Признаков: {len(combined_df.columns)}")
return combined_df
# =============================================================================
# ГЛАВНАЯ ФУНКЦИЯ
# =============================================================================
def main():
print_structure()
check_data()
if not DATA_ROOT.exists():
print("\n❌ Папка с данными не найдена!")
return
try:
extractor = MultiChannelRadiomicsExtractor()
print("\n" + "=" * 60)
print("📊 ИЗВЛЕЧЕНИЕ ПРИЗНАКОВ (ВСЕ КАНАЛЫ)")
print("=" * 60)
results = extractor.process_all_patients(
images_dir=PATH_IMAGES_TRAIN,
masks_dir=PATH_LABELS_TRAIN,
output_dir=PATH_FEATURES
)
if results:
print("\n" + "=" * 60)
print("📊 ОБЪЕДИНЕНИЕ РЕЗУЛЬТАТОВ")
print("=" * 60)
combined_df = merge_features(PATH_FEATURES)
if combined_df is not None:
visualize_features(combined_df, PATH_VISUALIZATION)
print("\n" + "=" * 70)
print("✨ ГОТОВО!")
print(f" 📂 Данные: {DATA_ROOT}")
print(f" 📂 Результаты: {OUTPUT_ROOT}")
print(f" ├── features/ {PATH_FEATURES}")
print(f" ├── visualization/{PATH_VISUALIZATION}")
print(f" └── logs/ {PATH_LOGS}")
print("=" * 70)
except Exception as e:
logger.error(f"❌ Ошибка: {e}")
import traceback
traceback.print_exc()
if __name__ == "__main__":
try:
main()
except KeyboardInterrupt:
print("\n⚠️ Прервано")
except Exception as e:
print(f"\n❌ Ошибка: {e}")
import traceback
traceback.print_exc()