Optimizing radiotherapy plans for cancer treatment with Tensor Networks解读

论文链接：https://arxiv.org/abs/2010.09552
本论文的主要内容（除去附录）

1.背景

放射治疗(后面称为放疗)：是利用电离辐射治疗实体肿瘤的技术之一
原理：释放到癌组织中的辐射剂量会损害肿瘤细胞的DNA，导致其死亡或减缓肿瘤的生长
用电离辐射治疗时最困难的挑战之一：

向目标肿瘤提供最佳剂量，同时保持周围健康组织中的辐射尽可能低

---

当今最常用的技术之一：强度调制放疗（IMRT）
IMRT的核心思想：对所施加的辐射进行调制，以在患者体内达到最佳剂量分布
缺点：由于该优化调制通过求解一个困难的数值优化问题获得的，所以有大量的优化参数和对最终辐射剂量分布的众多约束

---

在过去的几年中，已经发展了几种数值技术来应对这一挑战：

2015年，Nazareth和Spaans提出使用 D波退火机解决IMRT束流通量优化问题 2020年， Pakela等人提出了一种利用 模拟量子隧道效应的方法 这些新方法的应用有限： 1.缺乏可伸缩的量子硬件 2.缺乏将经典IMRT问题映射到量子硬件的稳健而清晰的策略 解决方案：将张量网络（TNs）应用于IMRT剂量优化问题

---

问题：如何使用TNs解决IMRT优化问题呢？

---

基本思路：

在本论文中将TNs应用于：
1.两个不同的简单模型(toy models)
2.一个更真实的模拟前列腺癌治疗的解剖场

本实验最终主要得到以下结果：
证明了TNs结果与其他经典技术(二次规划(QP)和模拟退火(SA))是相容( compatible )的

2.Methods

2.1 Radiotherapy optimization problem(放疗优化问题)

IMRT：在提高治疗质量方面，最具影响力的技术之一
在一些病例中，被关键OARs(危险区域：在病原细胞附近的会被影响的区域)包围的非常不规则的肿瘤需要高剂量，选择IMRT治疗因为它 辐射束(beam)非均匀强度分布
IMRT的原理：将每个beam的流量分割成一定数量的小pencil-beam，也称为𝑏𝑒𝑎𝑚𝑙𝑒𝑡(通过多叶准直器的移动传送)，从而对辐射beam进行调节

注：beam指的是一个角度射出的光，beamlet是对beam的分解，将其划分为一个一个的小格子

---

Challenge：找到每个beamlet的理想强度，从而优化患者体内的场分布和剂量分布之间的关系

---

该优化问题通常作为一个逆问题来解决：

虽然在过去几年中，提出了不同的方法：

但找到患者特定的方案仍然是一个挑战，因为永远都很难达到一个真正的理想剂量，同时它需要考虑许多因素，特别是随着时代的发展，要考虑的因素越来越多，从而增加了问题的复杂性。

IMRT的典型放疗程序

治疗计划中规定了tumor和OARs的预期剂量，目标是在向tumor提供适当剂量的同时保留健康器官

注：只优化穿过tumor的𝑏𝑒𝑎𝑚𝑙𝑒𝑡

患者要优化的区域可划分为三维有限尺寸单元𝑣𝑜𝑥𝑒𝑙𝑠

每个𝑏𝑒𝑎𝑚𝑙𝑒𝑡对每个𝑣𝑜𝑥𝑒𝑙的贡献取决于：

1.自身强度值 2.每个beam的几何和物理特性 3.患者的特征

---

所有这些信息都收集在𝑖𝑛𝑓𝑙𝑢𝑒𝑛𝑐𝑒 𝑚𝑎𝑡𝑟𝑖𝑥 𝐴中

𝑣𝑜𝑥𝑒𝑙 𝑖的总剂量可以表示为：

𝑁_𝐵: 𝑏𝑒𝑎𝑚𝑙𝑒𝑡的总数
𝑎_𝑖𝑗: influence matrix,A
注：𝑖𝑛𝑓𝑙𝑢𝑒𝑛𝑐𝑒 𝑚𝑎𝑡𝑟𝑖𝑥通常由临床环境中使用的治疗计划系统提供

---

通过迭代逆规划过程求解IMRT规划过程：
应用photon beams的强度朝着患者体内规定的剂量分布$D_i^{(P)}$进行优化
以beamlet权重x为优化变量，最小化患者身体中离散volume的交付剂量D(x)和期望剂量$D_i^{(P)}$之间的距离
通过以下二次成本函数来描述该过程：

r：覆盖所有volumes的索引
$v_r$属于对象r的voxels总数
${\gamma }_i$分配给第i个voxel的权重(以便在治疗期间优先考虑某些volumes)

2.2 Solving classical problems with TNs

2.2.1 Mapping the problem to a classical Hamiltonian

用一组$N_Q$ bits将每个beamlet 的权重$x_j$离散化:

$N_Q$ :位深度
$b_n^{(j)}$={0,1}
B:归一化常数
引入归一化常数B来设置beamlet权重的范围:

注：位数$N_Q$ ↑，离散化的分辨率↑
考虑经典哈密顿量 ${\mathrm{H}}_{IMRT}$ ，进一步允许每个自旋都是一个量子变量，用泡利矩阵${\sigma }_z$表示为$s_n^{(j)}$

构造了一个$N_B\times N_Q$维的多体哈密顿量
注：第一维在不同的beamlet上运行，第二维代表每个beamlet的离散化空间
得到伊辛型哈密顿量：


注：该公式的推导将在另一篇博客中展示https://blog.csdn.net/weixin_49883619/article/details/120685464

---

此公式相关信息：
1.上述方程的经典伊辛型哈密顿量的基态等价于解决最小化方程的初始问题
2.系统的基态由提供最低能量 $E_0$(使等式最小化的配置)的经典自旋组态给出
3.$E_0=0$对应于最佳的beamlet设置，从而在患者体内精确获得所需的剂量分布

---

在实践中，精确的理想剂量分布并不总是能够实现优化，可以通过伊辛型哈密顿量的所有可能自旋组态来找到初始优化问题的解
缺点：因为自旋构型的数量随着系统尺寸的增加呈指数增长，所以这种经典的搜索不可行

注：放疗所需的beamlet数量很容易达到几千个左右，因此代表它们所需的自旋数量实际上太高，无法在经典方案中有效解决问题

解决方案：基于TNs的方法来解决这个复杂的优化任务

2.2.2 Description of Tree Tensor Networks (TTN) for solving quantum many-body problems

在TTN处理后，得到的量子哈密顿量基态通常可能不是经典的，而是经典态的叠加

初始IMRT问题的所有解都是经典解，因此，一旦TTN算法完全收敛到量子基态，希望每个叠加的经典解
分别具有相同的基态能量$E_0$

生成的自旋进一步映射回每个voxel的二进制编码解决方案

结论：在大多数情况下，得到的TTN基态确实是一个经典态（$\chi =1$），或者在某些情况下，得到的态如预期的那样是经典解的叠加，每个经典解具有相同的能量$E_0$

注意：TTN算法对初始化过程非常敏感，对于每次运行，随机初始化几个TTN样本，从中开始优化，最后，通过比较得到的能量来验证最佳模拟

3.Results

3.1 Toy models

使用两种简单模型验证本文提出方法的准确性，和QP和SA两种优化方法对比，检验实验结果

3.1.1 Bipartite box

将三维盒子组成的模型分为两个区域，不同的区域使用不同的方式处理

红色数字表示每个部分的期望beamlet 权重
从两个相反的方向辐射立方体，这里表示为两个矩阵
每个beam的beamlets表示为$N_B$，都是可变的

将观察到它们的和等于每个二分的期望值，显示了解析解和数值解之间的完全一致。证明了算法的有效性和映射过程的正确性

3.1.2 Sphere

在0°、120°、180°和240°下使用四束光照射盒子，总数量$N_B$=64 beamlets（16 beamlets/beam）
𝛾 = 1.0

3.2. Prostate cancer treatment with TG119 IMRT phantom

该数据集包含多个分段结构，为了模拟前列腺癌的病例，我们选择了以下结构: 前列腺作为靶向肿瘤，膀胱和直肠作为OARs

在90°、270°处放置两个beams（33和31个beamlets）

4.Conclusions

注意：实验部分没有细讲，需要自己看论文，但是总的就是论证的本论文的方法与QP和SA方法具有较好的相容性