2024年认证杯SPSSPRO杯数学建模B题(第二阶段)神经外科手术的定位与导航全过程文档及程序
2024年认证杯SPSSPRO杯数学建模
B题 神经外科手术的定位与导航
原题再现：
人的大脑结构非常复杂，内部交织密布着神经和血管，所以在大脑内做手术具有非常高的精细和复杂程度。例如神经外科的肿瘤切除手术或血肿清除手术，通常需要将颅骨打开一个（或几个）圆形窗口，将病变部位暴露在术野中。但当病变部位较深时，就必须将上方的脑组织进行一定程度的牵拉、移位。这种牵拉和移位的幅度不能太大，而且具有损伤脑组织的风险。所以医生需要仔细设计手术方案，才能在尽量减小损伤的前提下完成手术目标。现在有CT或核磁共振等成像手段可以对颅内区域进行准确的透视和三维定位。但是在打开颅骨后，由于打开部位的颅骨约束作用消失，而且脑组织是柔软的，所以脑组织会产生一定程度的变形、移位，甚至在一定程度上膨出，以至于病变部位的位置发生了改变。如果想要在术中实时进行透视和定位，需要复杂的设备，术者也需经过特殊的训练，这都使得此类手术尚未完全普及。所以我们需要设计一个合理的数学模型，较为精细地预测手术部位在开颅后的位置改变情况。我们一定能够取得的数据是：
1. 术前的CT成像结果，这包含了颅腔的三维形状以及手术部位在三维空间中的位置（在建模时可自行寻找并使用典型的颅脑CT成像数据）；
2. 术前的颅内压；
3. 在颅骨上打开的窗口的大小以及位置；
4. 开颅后在颅骨窗口处可观察到的脑组织膨出高度（通常以mm计）;
5. 在术前可以测量得到的其他生理指标。
请你和你的团队建立合理的数学模型以解决下列问题。
第二阶段问题：
1. 假设患者的大脑中有一硬质占位性病变，我们通过术前的CT成像结果可以准确地知道它的位置。在术中可以测得颅内压。请尽可能精细地计算该占位性病变在开颅后的位移情况。
2. 假设患者的大脑中有一血肿区域（内容物为液态血液）。我们通过术前的CT成像结果可以准确地知道它的位置。在术中可以测得颅内压。请尽可能精细地计算该血肿在开颅后的位移和可能的形变情况。
整体求解过程概述(摘要)
神经外科手术，特别是涉及脑内深部病灶的手术，如肿瘤切除或血肿清除，面临着极高的技术挑战。这类手术要求在保护正常脑组织不受损伤的同时，精准定位并处理病灶。然而，开颅操作会导致颅内压变化，进而使脑组织发生难以预测的位移，影响术中病变位置的准确性，增加了手术难度和风险。为应对这一挑战，本文旨在通过综合考虑CT 图像、手术切口位置、脑组织力学特性、颅内压变化等因素，通过构建人工智能模型预测开颅后脑组织的位移情况。本文没有使用传统的有限元方法建模，有限元建模计算资源需求高、时间消耗长、非专业用户友好性极差，限制了其在临床环境中的广泛应用。而基于迁移学习的预训练-微调模型解决了这些问题，通过预先训练的模型共享参数后进行微调，极大的降低了运算资源以及数据需求量，对小样本问题友好；不需要丰富专业知识进行复杂物理建模分析，更好地服务于临床决策和手术规划。
我们分析了人工智能在医学图像领域的工作以及预测可行性，根据上述思想设计了通用医学影像预训练-微调模型，首先通过特征工程对105个特征进行处理，使用person相关系进行分析，通过向预训练模块输入颅脑等器官CT影像数据，使模型学习医学影像通用信息，共享预训练模块学习到的通用参数，再向模型输入少量肿瘤位移前后影像信息以及临床文本信息进行微调，使模型学习到更有针对性的细节信息。最后通过神经网络，输出目标转移到大脑某区域位置概率。
问题一，首先利用获取的CT图像获取占位性病变的精确尺寸、位置、形状以及与周围脑组织的关系，进行三维重建。对物理属性参数，如脑组织物理特性，粘弹性、颅内压、骨窗大小等特性进行分析，模拟硬质占位性病变位移范围，划定颅脑范围作为计算的输入参数。最后，通过医学影像预训练-微调模型预测位移概率。
问题二，血肿与问题一中的硬质占位性病变不同，可以简单理解为固体与液体的区别，我们通过分析液体动力学、固体力学以及脑组织的非线性弹性行为，分析了血肿在开颅手术后的位移和形变。根据脑组织的力学属性和血肿的流体力学属性设计了新的损失函数。通过修改微调信息，来预测血肿在开颅后的位移和形变。
最后，理想的神经外科手术导航系统应具备以下特点：提供高分辨率的实时影像引导，能够适应脑组织的动态变化，拥有精确的术中定位能力，并且便于外科医生操作，减少手术时间和并发症。随着计算机科学、图像处理技术、机器学习以及生物力学模型的不断进步，未来的神经外科手术有望实现更加个性化、精准的治疗方案，提高手术成功率并降低患者风险。
问题分析：
在分析本题的过程中我们注意到，整个脑部组织及其他人体特征的多样性。首先在题目中提到的数据，如表一所示。
针对数据问题：在Github的医学影像数据集集锦仓库1中，我们获得了具有各类模态以及各类器官约20个方向的80+个数据集。通过筛选选择了105个与本问题相关的参数，并且根据脑组的粘弹性等物理特性，开颅位置、大小以及膨出高度进行了数据模拟。
针对问题一：在各类器官影像数据量大，而相关肿瘤开颅位移影像数据量稀少的情况下，我们将颅内组织位置通过使用迁移学习的思路，利用人工智能的方法对硬质占位性病变在开颅后的位移情况进行预测。首先采取迁移学习中经典的预训练-微调模型，通过向预训练模块输入大量各种器官的CT影像数据，使预训练模块学习与CT影像相关的大量通用信息；其次，根据(Schirmer,Gieseetal.2019)的文章我们将全脑划分为10个不同的位置，我们将少量的包含肿瘤位移前后信息的三维数据以及文本信息数据，输入到模型中进行微调；最后，模型输出硬质占位性病变转移位置的概率，通过联合损失函数，留一交叉验证对测试数据进行验证，用ACC、F1分数、HausdorffDistance（豪斯多夫距离）对模型进行评估。
针对问题二：据脑组织的力学属性和血肿的流体力学属性，我们在问题一的基础上，仍然使用迁移学习进行预训练微调的方法，进行血肿位移和形变预测。在预训练模块不需要做任何改动，直接调用预训练模块的通用参数即可，对于微调模块，我们根据血肿的特性，通过分析液体动力学、固体力学以及脑组织的非线性弹性行为，设计了新的综合损失函数。最后通过模型预测血肿在开颅大致的位移和形变。根据已有的有限数据，我们可以通过计算机视觉模型来辅助判断出病变的类型。脑组织具有独特的生物力学特性，其中包括各向异性和粘弹性[1]，这些特性对脑部病变的诊断和研究至关重要。各向异性指的是脑组织在不同方向上的力学性能不同，而粘弹性则描述了脑组织在应力作用下的变形行为及其时间依赖性。通过对这些生物力学特性的分析和计算，我们可以粗略估计出肿瘤的位置偏移。
为了进一步提升预测的准确性和模型的泛化能力，我们采用迁移学习的方法。迁移学习是一种将已有知识应用到新的但相关问题上的技术。具体来说，我们使用目前已有关于出血性脑卒中的研究成果作为迁移学习的基础数据。这些研究成果包括大量的原始数据，这些数据经过充分的标注和分析，能够为我们的模型提供强大的先验知识基础。通过使用这些出血性脑卒中的原始数据，我们能够训练一个初始模型，该模型可以实现对血肿及水肿位置的预测。在此基础上，我们进一步将脑中病变及生理数据输入该模型，并对其进行微调。这一过程称为模型的微调（fine-tuning），旨在使模型更好地适应新的特定任务，从而提高其预测精度。通过这种方式，我们能够实现对病变位移的准确预测，为临床诊断和治疗提供有力的辅助工具。这一系列步骤展示了计算机视觉技术在医学影像分析中的巨大潜力，特别是在复杂脑部病变的诊断和治疗方面。
模型假设：
假设获取的全部数据是准确无误的；
假设脑组织的特性符合粘弹性这一材料属性的各项指标；
假设颅骨内各个关联组织之间无相对运动；
假设正常脑组织和病变组织均质地均匀；
假设忽略脑中的微观结构；
假设开颅过程是理想情况，一切操作均符合相关规范要求。
完整论文缩略图
全部论文及代码请见下方“ 只会建模 QQ名片” 点击QQ名片即可
部分程序代码：
import os
import torch
import torchvision as tv
import torchvision.transforms as transforms
import torch.nn as nn
import torch.optim as optim
import argparse
import skimage.data
import skimage.io
import skimage.transform
import numpy as np
import matplotlib.pyplot as plt
import torchvision.models as models
from PIL import Image
import cv2
#提取某一层网络特征图
class FeatureExtractor(nn.Module):
def __init__(self, submodule, extracted_layers):
super(FeatureExtractor, self).__init__()
self.submodule = submodule
self.extracted_layers = extracted_layers
def forward(self, x):
outputs = {}
for name, module in self.submodule._modules.items():
if "fc" in name:
x = x.view(x.size(0), -1)
x = module(x)
print(name)
if (self.extracted_layers is None) or (name in self.extracted_layers and 'fc' not in name):
outputs[name] = x
# print(outputs)
return outputs
def get_picture(pic_name, transform):
img = skimage.io.imread(pic_name)
img = skimage.transform.resize(img, (256, 256)) #读入图片时将图片resize成(256,256)的
img = np.asarray(img, dtype=np.float32)
return transform(img)
def make_dirs(path):
if os.path.exists(path) is False:
os.makedirs(path)
pic_dir = 'dataset/brain/train/1700.jpg'
transform = transforms.ToTensor()
img = get_picture(pic_dir, transform)
# 插入维度
img = img.unsqueeze(0)
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
img = img.to(device)
net = models.resnet101(pretrained=True).to(device)
dst = './feautures'
therd_size = 256
myexactor = FeatureExtractor(submodule=net, extracted_layers=None)
output = myexactor(img)
#output是dict
#dict_keys(['conv1', 'bn1', 'relu', 'maxpool', 'layer1', 'layer2', 'layer3', 'layer4', 'avgpool', 'fc'])
for idx,val in enumerate(output.items()):
k,v = val
features = v[0]
iter_range = features.shape[0]
for i in range(iter_range):
# plt.imshow(features.data.cpu().numpy()[i,:,:],cmap='jet')
if 'fc' in k: #不可视化fc层
continue
feature = features.data.cpu().numpy()
feature_img = feature[i, :, :]
feature_img = np.asarray(feature_img * 255, dtype=np.uint8)
dst_path = os.path.join(dst, str(idx)+'-'+k)
make_dirs(dst_path)
feature_img = cv2.applyColorMap(feature_img, cv2.COLORMAP_JET)
if feature_img.shape[0] < therd_size:
tmp_file = os.path.join(dst_path, str(i) + '_' + str(therd_size) + '.png')
tmp_img = feature_img.copy()
tmp_img = cv2.resize(tmp_img, (therd_size, therd_size), interpolation=cv2.INTER_NEAREST)
cv2.imwrite(tmp_file, tmp_img)
dst_file = os.path.join(dst_path, str(i) + '.png')
cv2.imwrite(dst_file, feature_img)