DRDMsig/Data_Engineering
0
1#coding:utf-8
2'''
3write by ygq
4create on 2025-09-04
5
6OASIS（Open Access Series of Imaging Studies） 是一个旨在向科研界免费提供脑部MRI数据的项目。本横断面（Cross-Sectional）数据集是其第一个版本，发布于2007年。
7OASIS-1 是横断面的，意味着它无法捕捉个体随时间的动态变化。对于研究疾病进展，后续的 OASIS-2 和 OASIS-3（纵向数据集）是更好的选择。
8
9OASIS-2，全称为 Longitudinal Multimodal Neuroimaging: Principal 150 Subjects，是 OASIS 项目发布的第二个核心数据集。顾名思义，它的核心特点是 纵向（Longitudinal）。
10
11核心目标：
12 研究正常衰老和阿尔茨海默病（AD）中的大脑结构随时间变化的模式。
13研究设计：
14 纵向研究。同一批受试者被多次扫描和评估，持续数年。
15样本量：
16 150 名年龄在 60 到 96 岁之间的受试者。
17人群组成：
18 所有 150 名受试者在首次扫描时都被诊断为认知正常（CDR = 0）。
19 在研究期间，部分受试者仍然保持认知正常，而另一部分则发展为痴呆（被临床诊断为可能患有阿尔茨海默病）。
20数据采集：
21 每名受试者进行了 至少 2 次 的访视会话（session），最多达到了 5 次。
22 每次访视之间的平均间隔时间约为 2.2 年，整个研究跨度最长超过 7 年。
23 每次访视都包括：3-4 次 T1 加权 MRI 扫描（在单次会话中完成，用于平均以提高信噪比）和详细的临床神经心理评估。
24数据内容：
25 与 OASIS-1 类似，包括原始 DICOM 图像、预处理后的 Analyze 格式图像，以及全面的临床认知评估数据。
26
27
28关键区别的详细解释
29横断面 vs. 纵向 (Cross-Sectional vs. Longitudinal)：
30 OASIS-1 像是在给一个城市的所有人在同一天拍一张照片。你可以比较年轻人和老年人、健康人和病人的区别，但看不到任何一个人是如何变老或生病的。
31 OASIS-2 像是挑选了150位健康的老年人，然后每年都给他们拍一张照片，持续好几年。这样你就能亲眼看到有些人如何慢慢地出现变化，最终生病。这对于理解疾病的过程至关重要。
32受试者群体的区别：
33 OASIS-1 包含了已经确诊的AD患者，非常适合训练一个模型来学习“AD大脑看起来是什么样”。
34 OASIS-2 的受试者起点都是健康的，这使得它成为研究疾病前驱期（即临床症状出现之前）的宝贵资源。你可以分析那些最终患病的人，在多年前其大脑是否就已经存在细微的、可检测的差异。
35数据分析方法的差异：
36 分析 OASIS-1 通常使用跨主体（cross-sectional）比较，例如比较AD组和正常对照组的平均海马体积。
37 分析 OASIS-2 则侧重于个体内部随时间的变化（within-subject change）。例如，为每个受试者计算其年化脑萎缩率，然后比较保持正常组和转化组之间的萎缩速率差异。这需要更复杂的纵向统计模型。
38
39
40 1. 人口统计学信息
41 性别（M/F）
42 用手习惯（Hand）（均为右利手）
43 年龄（Age）
44 教育程度（Educ）（1-5级）
45 社会经济地位（SES）
46
47 2. 临床评估
48 MMSE（简易精神状态检查）
49 CDR（临床痴呆评级：0=正常，0.5=非常轻微，1=轻度，2=中度）
50
51 3. 衍生解剖指标
52 eTIV：估计颅内容积
53 ASF：图谱缩放因子
54 nWBV：标准化全脑体积
55'''
56import os
57import glob,re
58import pandas as pd
59import SimpleITK as sitk
60import argparse
61import json
62from tqdm import tqdm
63from util import meta_data
64import util
65import numpy as np
66# from bert_helper import *
67
68import shutil
69##dataset_meta
70import warnings
71warnings.filterwarnings("ignore")
72meta_id_name='MRI ID'
73##性别（M/F）,用手习惯（Hand）（均为右利手）,年龄（Age）,教育程度（Educ）（1-5级）,社会经济地位（SES）,MMSE（简易精神状态检查）,CDR（临床痴呆评级：0=正常，0.5=非常轻微，1=轻度，2=中度）,eTIV：估计颅内容积,ASF：图谱缩放因子,nWBV：标准化全脑体积
74#META_COLUMN=['ID', 'M/F', 'Hand', 'Age', 'Educ', 'SES', 'MMSE', 'CDR', 'eTIV','nWBV', 'ASF', 'Delay']
75META_COLUMN=['Subject ID', 'Group', 'Visit', 'MR Delay', 'M/F', 'Hand','Age', 'EDUC', 'SES', 'MMSE', 'CDR', 'eTIV', 'nWBV', 'ASF']
76
77TASK_VALUE="segmentation"
78CLAMP_RANGE_CT = [-300,300]
79CLAMP_RANGE_MRI = None # MRI images threshold placeholder TBC...
80TARGET_VOXEL_SPACING=None
81
82# def find_metadata_files(path):
83# # for Cancer Image Archive (TCIA) dataset
84# search_pattern = os.path.join(path, '**', 'metadata.csv')
85# return glob.glob(search_pattern, recursive=True)
86
87def find_metadata_files(path):
88 # for Cancer Image Archive (TCIA) dataset
89 search_pattern = os.path.join(path, '*.csv')
90 return glob.glob(search_pattern, recursive=True)
91##added by yanguoqing on 20250527
92def find_image_dirs(path):
93 return os.listdir(path)
94
95##modify by yanguoqing on 20250527
96def load_dicom_images(folder_path):
97 reader = sitk.ImageSeriesReader()
98 dicom_names = reader.GetGDCMSeriesFileNames(folder_path)
99 reader.SetFileNames(dicom_names)
100 image = reader.Execute()
101 return dicom_names,image
102
103##added by yanguoqing on 20250527
104def load_dicom_tag(imgs):
105 reader = sitk.ImageFileReader()
106 # dicom_names = reader.GetGDCMSeriesFileNames(folder_path)
107 reader.SetFileName(imgs)
108 reader.ReadImageInformation() # 仅读取元信息，不加载像素数据
109 # metadata_keys = reader.GetMetaDataKeys()
110 tag=reader.Execute()
111 return tag
112
113def load_nrrd(fp):
114 return sitk.ReadImage(fp)
115
116##modify by yanguoqing on 20250904
117def load_raw_images(series_files):
118 '''
119 每个病例包含3到4种RAW的单次平扫MR
120 将多个分开的模态合并，构建第四个维度的数组，分别按照MPR-1,MPR-2...顺序存放
121 '''
122 reader = sitk.ImageSeriesReader()
123 reader.SetFileNames(series_files)
124 image = reader.Execute()
125 return image
126
127def save_nifti(image, output_path, folder_path):
128 # Set metadata in the NIfTI file's header
129 output_dirpath = os.path.dirname(output_path)
130 if not os.path.exists(output_dirpath):
131 print(f"Creating directory {output_dirpath}")
132 os.makedirs(output_dirpath)
133 # Set metadata in the NIfTI file's header
134 image.SetMetaData("FolderPath", folder_path)
135 sitk.WriteImage(image, output_path)
136
137##modify by yanguoqing on 20250527
138def convert_windows_to_linux_path(windows_path):
139 # Replace backslashes with forward slashes and remove the drive letter
140 # Some meta files have windows paths, but the data is stored on a linux server
141 linux_path = windows_path.replace('\\', '/')
142 if ':' in linux_path:
143 linux_path = linux_path.split(':', 1)[1]
144 return linux_path
145
146def main(target_path, output_dir):
147 pid_dirs=find_image_dirs(target_path)
148 failed_files = []
149 if not os.path.isdir(output_dir):
150 os.makedirs(output_dir)
151 json_output_path = os.path.join(output_dir, 'nifti_mappings.json')
152 failed_files_path = os.path.join(output_dir, 'failed_files.json')
153 meta = meta_data()
154
155 # Initialize the JSON file
156 if not os.path.exists(json_output_path):
157 with open(json_output_path, 'w') as json_file:
158 json.dump({}, json_file)
159 ##方便处理解析信息，转成csv文件
160 meta_file=os.path.join(os.path.dirname(os.path.realpath(__file__)),'oasis2_longitudinal_demographics.csv')
161 meta_file_ori=os.path.join(target_path,'oasis_longitudinal_demographics-8d83e569fa2e2d30 (1).xlsx')
162 if os.path.isfile(meta_file):
163 mf_flag=True
164 df_meta=pd.read_csv(meta_file,sep=',')
165 else:
166 mf_flag=False
167
168
169 if pid_dirs:
170 for pid_dir in tqdm(pid_dirs, desc="Processing pid dirs"):
171 if not os.path.isdir(os.path.join(target_path,pid_dir)):
172 continue
173
174 ##遍历所有目录下的病例数据
175 image_dirs=find_image_dirs(os.path.join(target_path,pid_dir))
176
177 for data_dir in tqdm(image_dirs, desc="Processing images files"):
178 ##data_dir即id
179 full_path=os.path.join(target_path,pid_dir,data_dir)
180
181 modality="MRI"
182 study='OASIS_2'##Dataset_name
183 CIA_other_info = {'metadata_file':''}
184 CIA_other_info['split'] = "train"
185 CIA_other_info['metadata_file']=meta_file_ori
186 data_info_row=df_meta[df_meta[meta_id_name]==data_dir]
187
188 if data_info_row.shape[0]>0:
189 data_info_row=data_info_row.reset_index()
190 #print(data_info_row[meta_id_name])
191 for keyname in META_COLUMN[:]:
192 CIA_other_info[keyname]=str(data_info_row[keyname][0])
193
194 CIA_other_info['Image_id']=data_dir
195
196
197 else:
198 meta_image_id=data_dir
199 for keyname in META_COLUMN[1:]:
200 CIA_other_info[keyname]=''
201
202
203
204 try:
205 ##读取原始的RAW目录下多次单扫img
206 #\RAW\OAS1_0001_MR1_mpr-1_anon.img
207 series_files=glob.glob("%s/RAW/mpr-*.img"%(full_path))
208 series_files.sort()
209
210 if len(series_files)>0:
211 ##存在有效的MRI影像数据进行后续处理
212 sitk_img_original=load_raw_images(series_files)
213 submodality=[re.search(r"mpr-\d{1}",os.path.basename(fp)).group(0) for fp in series_files]
214 sub_modality_dict={}
215 for idx,value in enumerate(submodality):
216 sub_modality_dict[idx]=value
217
218 meta.add_keyvalue('Sub_modality',sub_modality_dict)
219
220 else:
221 print("病例数据%s为空"%data_dir)
222 continue
223
224 original_spacing = list(sitk_img_original.GetSpacing())
225 original_size = list(sitk_img_original.GetSize())
226 print(original_spacing)
227 is_4d_image = sitk_img_original.GetDimension() == 4
228
229 frame_flag=False
230 # --- Resampling Logic (Revised for 4D) ---
231 if is_4d_image:
232
233 # Always process 4D images channel-wise for resampling
234 # logging.info(f" Processing 4D image channel-wise: {original_img_full_path}") # Keep log for errors only
235 channels = []
236 num_channels = original_size[3] if len(original_size) == 4 and sitk_img_original.GetDimension() == 4 else 1
237 channel_target_spacing = TARGET_VOXEL_SPACING if TARGET_VOXEL_SPACING else original_spacing[:3] # Use 3D spacing
238
239
240 for i in range(num_channels):
241 extractor = sitk.ExtractImageFilter()
242 current_3d_channel_size = original_size[:3]
243
244 if sitk_img_original.GetDimension() == 4:
245 extractor.SetSize([current_3d_channel_size[0], current_3d_channel_size[1], current_3d_channel_size[2], 0])
246 extractor.SetIndex([0,0,0,i])
247 channel_3d_img = extractor.Execute(sitk_img_original)
248 else:
249 channel_3d_img = sitk_img_original
250 if i > 0: break
251
252 channel_resampler = util.get_unisize_resampler(
253 channel_3d_img, 'linear',
254 spacing=channel_target_spacing, size=current_3d_channel_size
255 )
256 if channel_resampler:
257 channels.append(channel_resampler.Execute(channel_3d_img))
258 else:
259 channels.append(channel_3d_img)
260
261 if channels:
262 if len(channels) > 1: # Only join if there are multiple channels
263 sitk_img_processed = sitk.JoinSeriesImageFilter().Execute(channels)
264 ##aded by yanguoqing on 2025-08-11
265 frame_flag=True
266 # imgDict={}
267 # for kf_idx in range(num_channels):
268 # imgDict[str(kf_idx)]='none'
269 # if str(meta_ed):imgDict[str(meta_ed)]='ed'
270 # if str(meta_es):imgDict[str(meta_es)]='es'
271 # meta.add_keyvalue('ImgDict',imgDict)
272 elif len(channels) == 1: # If only one channel resulted (e.g. original was 3D misidentified as 4D by tensorImageSize)
273 sitk_img_processed = channels[0]
274 elif TARGET_VOXEL_SPACING: # 3D image with target spacing
275 img_resampler_obj = util.get_unisize_resampler(sitk_img_original, 'linear',
276 spacing=TARGET_VOXEL_SPACING, size=original_size)
277 if img_resampler_obj: sitk_img_processed = img_resampler_obj.Execute(sitk_img_original)
278 else: # 3D image, no TARGET_VOXEL_SPACING
279 img_resampler_obj = util.get_unisize_resampler(sitk_img_original, 'linear',
280 spacing=original_spacing, size=original_size)
281 if img_resampler_obj: sitk_img_processed = img_resampler_obj.Execute(sitk_img_original)
282
283
284
285
286
287 original_spacing = list(sitk_img_original.GetSpacing())
288 original_size = list(sitk_img_original.GetSize())
289 size_processed = list(sitk_img_processed.GetSize())
290 print('size_processed',size_processed,original_size)
291
292
293 meta.add_keyvalue('Spacing_mm',min(original_spacing[:3]))
294 meta.add_keyvalue('OriImg_path',",".join(series_files))
295 meta.add_keyvalue('Size',size_processed) # 这里用处理后的size -- YH Jachin
296 meta.add_keyvalue('Modality',modality)
297 meta.add_keyvalue('Dataset_name',study)
298 meta.add_keyvalue('ROI','head')
299
300
301
302 output_image_file = os.path.join(output_dir,data_dir, f"{data_dir}.nii.gz")
303 # output_path=convert_windows_to_linux_path(output_path)
304 ##
305 save_nifti(sitk_img_processed, output_image_file, full_path)
306 print(f"Saved NIfTI file to {output_image_file}")
307 ##Label processing
308
309
310
311 except Exception as e:
312 print(e)
313 failed_files.append(data_dir)
314 print(f"Failed to load OASIS images from {data_dir}")
315 continue
316
317
318
319 meta.add_extra_keyvalue('Metadata',CIA_other_info)
320
321
322 # Write the mapping to the JSON file on the fly
323 with open(json_output_path, 'r+') as json_file:
324 existing_mappings = json.load(json_file)
325 existing_mappings[output_image_file] = meta.get_meta_data()
326 json_file.seek(0)
327 # print(existing_mappings)
328 json.dump(existing_mappings, json_file, indent=4)
329 json_file.truncate()
330 # else:
331 # print("No metadata.csv files found.")
332
333 with open(failed_files_path, "w") as json_file:
334 json.dump(failed_files, json_file)
335
336 print(f"The list has been written to {failed_files_path}")
337 print(f"Saved NIfTI mappings to {json_output_path}")
338
339if __name__ == "__main__":
340 parser = argparse.ArgumentParser(description="Process DICOM files and save as NIfTI.")
341 parser.add_argument("--target_path", type=str, help="Path to the target directory containing metadata files.", default="/home/data/Github/data/data_gen_def/DATASETS/OASIS/OASIS_2/OAS2_RAW//")
342 parser.add_argument("--output_dir", type=str, help="Directory to save the NIfTI files.", default="/home/data/Github/data/data_gen_def/DATASETS_processed/OASIS/OASIS_2/RAW_V2")
343 args = parser.parse_args()
344 print(args.target_path, args.output_dir)
345 main(args.target_path, args.output_dir)