青光眼是全球排名首位的不可逆致盲眼病，定期的视野检查是青光眼诊断和治疗过程中的必要监测手段，提前预测患者未来视野将有利于临床医生对病情进展进行及时干预。为了联合利用患者过去视野检查结果中的时间和空间特征，以提高视野预测效果，本文采用卷积长短期记忆（ConvLSTM）网络构建预测模型，并使用来自华盛顿大学汉弗瑞视野分析仪的视野测试数据集（UWHVF）的数据，对ConvLSTM模型与其他方法进行预测性能验证与比较。研究结果显示，相较于传统方法，ConvLSTM模型具有更高的预测精度；同时，探究视野序列长度与预测性能的变化关系发现，当采用过去1.5～6.0年内的3次视野结果预测时，ConvLSTM模型的预测性能更好，预测结果的平均绝对误差为2.255 dB，均方根误差为3.457 dB，决定系数为0.960。实验结果表明，本文所提方法仅使用既往视野检测结果，即实现了较准确的未来0.5～2.0年内的视野预测，因此该方法有望用于辅助临床医生对视野进展进行评估并治疗。
引用本文： 王握, 郑秀娟, 吕智清, 李妮, 陈俊. 基于时空特征学习的视野预测研究. 生物医学工程学杂志, 2024, 41(5): 1003-1011. doi: 10.7507/1001-5515.202310072 复制
0 引言
青光眼是一种进展缓慢、损害视神经及视觉通路并最终导致视觉功能损害的疾病，是全球首位的不可逆性致盲病因[1]。目前中国青光眼患者超过2 100万，约占全世界青光眼患者的四分之一[2-3]。青光眼的特征是视网膜神经节细胞的丢失、视网膜神经纤维层（retinal nerve fiber layer，RNFL）变薄，并出现特征性视神经萎缩，而这些改变在功能上表现为特征性视野（visual field，VF）缺损[4]。视野改变是诊断青光眼的“金标准”，定期的视野随访检查是青光眼诊断和监测病情进展的重要手段，且该监测往往是终身持续的。由于青光眼患者的视野变化波动比正常人更严重[5]，视野恶化速度因个体而异，其中稳定型患者恶化速度通常小于0.5 dB/a，而快速进展型患者则超过1.5 dB/a[6]。并且随着时间的推移，患者的视野进展模式存在明显不同[7]，再加上患者常常无法保证定期检查，因此提前预测未来视野变化对青光眼诊断和治疗具有重要意义。特别是预测未来0.5～2.0年内视野变化，对医生尽早发现快速进展型患者，及时调整治疗方案，具有重要临床价值[8]。
目前，视野预测还主要聚焦于对视野进展评估研究，以期直接辅助临床医生提前做出治疗决策。视野进展评估是一种分类任务，但是对于青光眼的视野进展的定义目前缺乏一种公认的“金标准”。常用的视野进展定义大致有4类：第一类是基于专家经验的定义，即有经验的临床医生判断视野进展[9-10]；第二类是基于随机临床试验中采用的分级评分系统的定义，即根据视野评分进行进展分类[11-12]；第三类是基于事件的定义，以视野随访与基线检测的对比结果作为进展标准[9, 13]；第四类是基于视野参数变化趋势的定义，一般通过参数变化斜率及其统计学指标P值来判断视野是否进展[14-15]。这些视野进展的定义方法各有不同，多种定义方法的比较研究结果表明：不同定义得到的视野进展的分类结果存在较大差异[16-18]。
由于各类视野进展定义方法本质上都是基于视野各点敏感度数值的计算，如果能直接对未来视野的逐点敏感度进行预测，预测结果将不受临床医生选择视野进展定义的限制。目前对视野逐点敏感度进行预测的研究工作较少。Kamalipour等[19]采用一种卷积神经网络（convolutional neural network，CNN）模型，利用谱域光学相干断层扫描，将得到的患者RNFL厚度值作为输入，预测视野的逐点敏感度。Taketani等[20]采用了多种传统回归模型，如普通最小二乘线性回归（linear regression，LR）、二次回归和指数回归等，对视野逐点敏感度进行预测，最终发现要达到精确预测未来视野结果，所需的最少视野检查结果数目大约为10。Wen等[21]开发了一种CNN模型5层级联网络（CascadeNet5），该模型仅使用单个视野结果作为输入，即可预测青光眼患眼未来5.5年内的视野逐点敏感度。但是这种模型仅利用单个视野结果中的空间信息，而未充分利用患者过去随访结果中的时间信息。与Wen等[21]的工作不同，Park等[22]的工作更关注于视野随访结果中的时序信息，他们建立了一种递归神经网络（recurrent neural network，RNN）模型，该模型以一系列5个连续的视野结果作为输入，来预测第6个视野结果，该研究预测的是未来视野的总偏差（total deviation，TD）值，而不是逐点敏感度。
对于视野逐点敏感度预测而言，仅基于时间或空间单因素的预测方法可能造成过去视野检查结果中有效信息的浪费，而卷积长短期记忆（convolutional long-short term memory，ConvLSTM）网络是基于长短期记忆（long-short term memory，LSTM）网络的一种改进模型，它可以同时利用数据中的时间信息以及空间信息[23]。已有多项研究证明，在时空序列预测问题上，ConvLSTM模型可表现出相当好的性能[24-26]。因此，为实现更高精度的未来视野预测，本研究采用ConvLSTM网络作为预测模型。本文仅以患者过去视野的随访结果序列作为输入，而无需其它生理结构信息，并将未来0.5～2.0年内视野的逐点敏感度预测作为输出；另外，以样本中包含视野结果的数目作为视野序列长度，对不同长度视野序列对预测结果产生的影响进行了探究，以期找到预测效果最好的视野序列长度。
1 方法
1.1 ConvLSTM网络
本研究旨在联合利用患者过去视野序列的时序特征以及单个视野结果中各点敏感度的空间分布特征，实现对未来视野随访的预测。传统LSTM在输入到状态、状态到状态的转换中采用了全连接结构，无法有效提取空间信息[27]。然而，ConvLSTM将LSTM与卷积操作相结合，将LSTM中的转换结构改进为卷积结构，使得网络具备了对空间结构特征的提取能力，这种方法更适用于兼具长短程时间特性和空间结构特性的视野逐点敏感度预测任务。而且，ConvLSTM在每个时间步都使用相同的卷积核进行卷积操作，实现了参数共享和稀疏连接，相较于全连接结构其参数规模更小，运算效率更高。
ConvLSTM的结构与传统的LSTM基本相似，同样包括遗忘门（F）、输入门（I）、细胞状态（C）和输出门（O）。ConvLSTM的细胞结构示意如图1所示，各部分之间的关系如下：
图1
ConvLSTM细胞结构
Figure1.
ConvLSTM cell structures
（1）遗忘门：它决定了细胞状态中有多少部分被遗忘或保持，其计算过程如式（1）所示：
 |
其中，
表示卷积；σ(∙)表示S型生长曲线（sigmoid）映射函数；Xt、Ft分别代表第t时间步的输入和遗忘门输出；Ht−1代表第t − 1时间步的隐藏状态输出；Kfx、Kfh分别代表遗忘门在隐藏状态和输入上的卷积核参数，bf为遗忘门对应的偏置参数。
（2）输入门：它决定了当前时间步的新信息有多少部分用于更新细胞状态，其计算过程如式（2）所示：
 |
其中，It代表第t时间步的输入门输出，Kix、Kih分别代表输入门在隐藏状态和输入上的卷积核参数，bi为输入门对应的偏置参数。
（3）细胞状态：它结合遗忘门输出和输入门输出来更新自身，其计算过程如式（3）所示：
 |
其中，
表示哈达玛（Hadamard）积，tanh(∙)为双曲正切函数；Ct − 1、Ct分别代表第t − 1和第t时间步的细胞状态；Kcx、Kch分别代表细胞状态在隐藏状态和输入上的卷积核参数，bc为细胞状态对应的偏置参数。
（4）输出门：它决定细胞状态的哪一部分将作为当前时间步的隐藏状态输出，其计算过程如式（4）～式（5）所示：
 |
 |
其中，Ot、Ht分别代表第t时间步的输出门输出和隐藏状态输出；Kox、Koh分别代表输出门在隐藏状态和输入上的卷积核参数，bo为输出门对应的偏置参数。
神经网络在广泛的具有挑战性的预测问题上取得成功的原因通常归结于其深度[28]。因此，本研究采用堆叠式ConvLSTM模块，以增加网络深度，从而获得视野序列在不同时间尺度上的特征表示。本研究所采用模型的网络结构包括1个堆叠式ConvLSTM模块和1个卷积层。其中，堆叠式ConvLSTM模块仅包含2个ConvLSTM单元，并且为避免梯度消失，在2个ConvLSTM单元之间添加了1个批标准化（batch normalization，BN）层，对前1个ConvLSTM单元所提取特征的分布进行调整。该模型结构图如图2所示，其中，nVFs表示长度为n的视野序列输入，VFn+1表示模型的预测输出为第n+1个视野结果。同时，n也与ConvLSTM的时间步数对应。由于以时序图表示的单个视野结果尺寸较小，选择较小的卷积核尺寸更适合用于提取视野特征，故设置2个ConvLSTM单元的卷积核大小为3×3，数目为32。最后的卷积层亦选择卷积核大小为3×3，数目为1，以减少通道数目，进而输出视野预测结果。
图2
ConvLSTM网络模型
Figure2.
ConvLSTM network model
1.2 对比模型
为了验证ConvLSTM模型在视野逐点敏感度预测方面的优势，本研究首先构建了4个基线模型，以提供基线预测性能。
本研究选取了2种高斯分布模型（Gaussian distribution model，GDM）和2种LR方法作为基线模型用于后续模型性能比对。其中，第一种GDM模型（记为GDM1）是基于Heijl等[29]工作中的视野进展速率的平均值和标准差。具体而言，本研究使用早期明显青光眼的视野敏感度进展速率数据，即均值为−0.36 dB/a，标准差为0.60 dB/a，将视野序列中逐视野预测结果取平均作为未来视野的最终预测。第二种GDM模型（记为GDM2）则是基于由训练集中统计得到的视野进展速率的平均值和标准差，之后按与GDM1相同步骤进行视野预测[21]。第一种LR模型（记为LR1）针对于患者群体，对视野各点敏感度进行LR操作，得到一组回归系数，以之预测视野。第二种LR模型（记为LR2）则针对于患者个体的逐点敏感度，以患者之前的随访数据预测其未来视野[22]。然后，本研究以ConvLSTM模型结构为基础，构建CNN模型和LSTM模型作为对比模型。CNN模型以单个视野结果作为输入，仅利用视野逐点敏感度的空间位置信息对未来视野进行预测。而LSTM模型则以视野序列为输入，利用序列中的时序信息预测未来视野。CNN模型的结构与ConvLSTM模型基本一致，它由3个卷积层组成，前2个卷积层的卷积核大小、卷积核数目等参数设置与ConvLSTM单元的参数设置基本相同，最后一层卷积层与ConvLSTM模型的最后一层卷积层完全相同。而LSTM模型中同样仅包含2个LSTM单元，并且它以全连接层代替最后的卷积层。
2 实验
本研究使用的开源数据集为来自华盛顿大学汉弗瑞视野分析仪的视野测试数据集（a dataset of perimetry tests from the Humphrey field analyzer at the University of Washington，UWHVF）[30]，在此基础上开展研究工作，以ConvLSTM模型与其他几种对比方法进行预测性能验证与比较。
2.1 数据集
UWHVF数据集中所有患者都接受了视野测量与评估软件汉弗瑞视野分析仪24-2程序模式（Humphrey Field Analyzer，24-2 Program）（Carl Zeiss Meditec Inc.，美国）的视野分析，使用瑞典交互式阈值算法（标准或快速）或全阈值策略执行。数据集包含了3 871位患者的7 428只眼睛的视野数据，共计28 943个视野，所有视野数据均采集自1998年—2018年。每只眼睛的每个视野在时间上进行了对齐，第1个视野被视为每只眼睛的基线视野。对所有纳入数据集的患者进行统计，以“中值[四分位范围]”形式显示：年龄为64[54, 73] 岁；随访的时间为2.49[1.11, 5.03] 年；基线视野的平均总偏差（mean total deviation，MTD）为−4.51[−8.01, −2.65] dB；模式标准差（pattern standard deviation，PSD）为2.41[1.70, 5.34] dB。UWHVF数据集中患者视野示例如图3所示，逐点敏感度精确到百分位，灰度浅深与敏感度大小对应[30]。目前在UWHVF数据集上开展的视野相关的研究较少，尤其在视野敏感度预测方面仍无相关工作发表。
图3
青光眼患者视野与正常视野对比示例
Figure3.
An example of patient's visual field with glaucoma compared to normal person's visual field
2.2 数据预处理
首先，在患者水平上，本研究将UWHVF数据集划分为两部分：80%作为训练集，用于训练和验证阶段；20%作为测试集。两个数据集中患者的统计信息如表1所示（表格仅统计信息完整的数据），可以发现，除性别信息部分缺失导致的差异以外，两个数据集中的数据分布基本相近。为了探究不同长度的视野序列对未来视野的预测性能，本研究将两个集合当中的视野按照时间间距与序列长度分别做了筛选，以得到不同长度的视野序列数据集。具体而言，筛选的标准为：在同一视野序列样本中，相邻视野检查的时间间隔必须在0.5～2.0年内。此区间的选定是考虑到青光眼视野进展较缓慢，时间间距小于0.5年视野变化可能不明显，不利于捕捉其变化规律，而当时间间距大于2.0年，则可认为是视野检测失访，不宜纳入时序研究。此外，为了尽可能多地获得样本，允许视野序列部分重复。例如，某患者的一只眼睛有7个视野测试结果，分别在33.75、34.21、35.00、36.52、37.43、38.65、41.02岁，现可将其分为2个序列长度为4的视野样本：（34.21、35.00、36.52、37.43岁）和（35.00、36.52、37.43、38.65岁），或者1个序列长度为5的样本：（34.21、35.00、36.52、37.43、38.65岁）。并且，部分样本的视野序列中包含全零视野，说明该患者已经造成不可逆致盲，无需纳入研究而予以剔除。
表1
训练集与测试集统计信息
Table1.
Statistics information of training and testing sets
最终，训练集和测试集中样本数目统计如表2所示。由于当序列长度大于7时，数据集中符合筛选条件的样本数目过少，故研究仅针对长度为3～7的视野序列样本（3VFs～7VFs）。对于序列长度为n+1的数据集，实验将以前n个视野结果作为输入序列预测第n+1个视野，记为nVFs→VFn+1，其中n=2, 3 …, 6。
表2
不同长度视野序列数据数目统计
Table2.
Statistics on the number of visual field series data of different lengths
在视野测量与评估软件汉弗瑞视野分析仪24-2程序模式（Carl Zeiss Meditec Inc.，美国）测量下，共获得54个视野测试点。剔除2个生理盲点后，仅使用剩下的52个测试点。每个视野被表示为8×9的矩阵，以保持其空间关系。对所有输入特征采用零均值化方法，将其标准化成均值为0，标准差为1。
2.3 实验参数与评价指标
研究实验平台采用64位Windows10操作系统（Microsoft Corp.，美国），编程语言为Python 3.6.0（Python Software Foundation，美国），而Python中的深度学习算法在TensorFlow 1.7.0（Google Inc.，美国）上使用开源神经网络平台Keras 2.2.0（Francois Chollet et al.，美国）实现。训练深度学习模型的迭代次数均设置为25，批大小为32，选择自适应矩估计（adaptive moment estimation，Adam）作为优化器，损失函数均采用均方误差（mean squared error，MSE）损失。采用5折交叉验证方式对这些模型进行训练，选用交叉验证中表现最优的模型参与后续在测试集上的预测性能对比。
预测模型的评价标准采用平均绝对误差（mean absolute error，MAE）、均方根误差（root mean squared error，RMSE）和决定系数R2，以避免不同评价标准本身存在的局限性。它们对应的计算如式(6)～式(8)所示：
 |
 |
 |
其中，
为第i个视野点的真实敏感度，
为模型预测第i个视野点的敏感度，
为真实敏感度均值，m为预测逐点敏感度的数目。
3 实验结果及分析
3.1 模型对比结果
为了分析ConvLSTM模型在视野预测任务中的表现，使用相同的数据集训练6种对比模型以及来自文献[21]和文献[22]的2种算法模型。其中，文献[22]提出的算法模型是一种RNN模型，因此无需修改模型结构即可接受不同长度视野序列输入。随后，利用各模型在测试阶段的预测结果计算不同的评价指标。结果统计如表3所示，最优值以加粗形式表示。
表3
不同模型预测性能对比
Table3.
Comparison of the prediction performance of different models
通过表3中各个模型的对比结果可以发现：① 从3项评价指标来看，5种深度学习方法的性能显著优于各基线模型，表明它们对视野敏感度变化规律的刻画能力更强，其中，性能最好的是ConvLSTM模型。与仅考虑视野空间或时间单因素信息的CNN和LSTM模型相比，ConvLSTM模型通过联合考虑两种信息，有效提高了视野预测的精度；② 两种GDM模型的性能基本优于两种LR模型，说明使用线性方式难以准确刻画视野敏感度的变化，因此在表现上不如基于数据分布的方法；③ GDM1模型相较于GDM2模型表现更佳，这表明仅依靠训练集得到的统计特征不足以全面描述患者群体特征，因此GDM1在测试集上的预测效果更好；④ LR2模型的表现优于LR1模型，后者的决定系数皆小于0.1，这表明不同患者个体的视野变化差异较大，基于患者整体的LR方法并不适用于视野预测。⑤ 来自文献[21]的CascadeNet5模型与构建的CNN模型相比性能更好，这可能是因为前者的模型结构更加复杂，并且在输入中加入了年龄信息。然而，ConvLSTM模型参数规模远小于CascadeNet5模型，却取得了更好的预测结果，说明时空特征的考虑比增加模型复杂度更加重要。⑥ 来自文献[22]的RNN模型在5种深度学习模型中表现最差，这可能是由于它最初设计用于预测视野的总偏差值而非逐点敏感度。此外，在输入的视野序列数据中人为加入了空白数据，增加了噪声，进一步提升了模型的学习难度。
对于不同缺损程度的视野序列，各深度学习模型的预测性能可能会受到不同程度的影响。本研究将5种深度学习模型应用于不同缺损程度的视野序列进行预测，如图4所示为代表性预测示例，其中包含预测输出值及其误差绝对值两种图像，缺损程度由平均总偏差确定。从图4中示例可以看出，ConvLSTM模型受到视野缺损程度影响相对其他模型更小，主要的预测偏差出现在敏感度数值“分界线”上。进一步地，对ConvLSTM模型在各测试集上的预测结果按照平均总偏差值分段进行统计，采用RMSE指标，并以箱线图方式表达，如图5所示。从图5中的对比结果可以看出，模型对于平均总偏差值较大的测试样本，预测结果基本较好且相对稳定。这表明，对于视野缺损严重的患者，模型能更好地预测其视野敏感度变化，有助于医生对病情发展进行及时干预。此外，由于各测试集中包含的样本数目差距较大，箱线图中“异常值”数目也与之相关，3VFs测试集上的“异常值”最多，7VFs测试集上最少。7VFs测试集中平均总偏差值为“>18 dB”的样本RMSE分布较分散，也是由于样本数目较少所导致。
图4
各模型对不同缺损程度视野的预测代表性示例
Figure4.
Representative examples of predicted visual fields with different defect degrees by different models
图5
ConvLSTM模型在不同序列长度测试集上的预测结果统计（按平均总偏差分段）
Figure5.
Statistics of prediction results of ConvLSTM model on different series length testing sets (segmented by mean total deviation)
3.2 不同序列长度结果
对于时序深度学习模型而言，其输入时间序列的时间步数可能影响最终的预测性能。为了探究不同长度的视野序列输入的影响，本研究对ConvLSTM模型在各测试集上的预测性能进行了比较分析。
由表3中的统计结果可以看出，视野序列长度与模型预测性能之间存在关系。对于ConvLSTM模型而言，结果显示视野序列长度并不是越长或越短预测性能越好，而是长度为3（即3VFs→VF4）的视野序列预测效果最佳。将LSTM和ConvLSTM两种时序深度学习模型的评价指标差值以折线图表示，如图6所示，其中Δ表示差值。可以看到，随着输入视野序列长度变化，MAE、RMSE和R2等指标的差值不断增大，说明二者性能差距在逐步拉开，而且LSTM模型预测指标原值是一直下降的，这表明视野序列中的时序信息的作用在减弱，而空间信息的作用则有所增强。原因在于，随着视野序列长度的增加，视野检测时间跨度增大，带来的时间冗余信息随之不断增多，因此空间信息在视野预测中所占比重上升。在6VFs→VF7的预测中，CNN模型和文献[21]模型预测结果的MAE值反而小于ConvLSTM模型，而其余指标均与后者相差不大，这说明在时间冗余信息的影响下，空间信息所起的作用在逐步凸显。
图6
LSTM和ConvLSTM模型的预测性能差距
Figure6.
Predictive performance gap between LSTM and ConvLSTM models
4 结语
目前，针对视野逐点敏感度预测方面开展的工作还很有限，而现有研究方法仍存在诸多不足，如预测精度低、视野序列中的时空信息利用不充分等。本研究采用ConvLSTM模型对未来0.5～2.0年内视野的逐点敏感度进行预测，以青光眼患者既往视野随访结果作为输入序列，而无需利用其他眼科特检检查或眼部解剖特征，实现了较高精度的预测。主要研究结论如下：① 相较于作为基线的传统模型，深度学习模型在刻画视野逐点敏感度变化方面的能力更强；② ConvLSTM模型充分考虑了视野敏感度随时间变化规律以及各位置点之间的空间关系，因此其预测精度高于CNN模型和LSTM模型；③ 输入预测模型的视野序列长度并非越长越好，时间序列过长可能带来更多冗余信息，实验发现适宜于ConvLSTM模型的输入序列长度为3，即采用过去1.5～6.0年内的3次视野检查结果预测效果更好。
当然，本研究在数据集选用方面存在一定局限性，由于目前开源的视野数据少，研究仅采用了UWHVF数据集，而缺少了其他人群的视野数据。同时，本研究工作是统计学意义上的改进，对于其能否在临床上进行应用，需要结合临床进一步评估，这也将是后续研究需要开展的工作。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：王握负责模型算法的设计和论文的撰写；郑秀娟负责指导方法理论和算法设计，以及指导论文的撰写；吕智清负责资料收集和论文修改；李妮和陈俊负责指导医学基础和分析方法。
0 引言
青光眼是一种进展缓慢、损害视神经及视觉通路并最终导致视觉功能损害的疾病，是全球首位的不可逆性致盲病因[1]。目前中国青光眼患者超过2 100万，约占全世界青光眼患者的四分之一[2-3]。青光眼的特征是视网膜神经节细胞的丢失、视网膜神经纤维层（retinal nerve fiber layer，RNFL）变薄，并出现特征性视神经萎缩，而这些改变在功能上表现为特征性视野（visual field，VF）缺损[4]。视野改变是诊断青光眼的“金标准”，定期的视野随访检查是青光眼诊断和监测病情进展的重要手段，且该监测往往是终身持续的。由于青光眼患者的视野变化波动比正常人更严重[5]，视野恶化速度因个体而异，其中稳定型患者恶化速度通常小于0.5 dB/a，而快速进展型患者则超过1.5 dB/a[6]。并且随着时间的推移，患者的视野进展模式存在明显不同[7]，再加上患者常常无法保证定期检查，因此提前预测未来视野变化对青光眼诊断和治疗具有重要意义。特别是预测未来0.5～2.0年内视野变化，对医生尽早发现快速进展型患者，及时调整治疗方案，具有重要临床价值[8]。
目前，视野预测还主要聚焦于对视野进展评估研究，以期直接辅助临床医生提前做出治疗决策。视野进展评估是一种分类任务，但是对于青光眼的视野进展的定义目前缺乏一种公认的“金标准”。常用的视野进展定义大致有4类：第一类是基于专家经验的定义，即有经验的临床医生判断视野进展[9-10]；第二类是基于随机临床试验中采用的分级评分系统的定义，即根据视野评分进行进展分类[11-12]；第三类是基于事件的定义，以视野随访与基线检测的对比结果作为进展标准[9, 13]；第四类是基于视野参数变化趋势的定义，一般通过参数变化斜率及其统计学指标P值来判断视野是否进展[14-15]。这些视野进展的定义方法各有不同，多种定义方法的比较研究结果表明：不同定义得到的视野进展的分类结果存在较大差异[16-18]。
由于各类视野进展定义方法本质上都是基于视野各点敏感度数值的计算，如果能直接对未来视野的逐点敏感度进行预测，预测结果将不受临床医生选择视野进展定义的限制。目前对视野逐点敏感度进行预测的研究工作较少。Kamalipour等[19]采用一种卷积神经网络（convolutional neural network，CNN）模型，利用谱域光学相干断层扫描，将得到的患者RNFL厚度值作为输入，预测视野的逐点敏感度。Taketani等[20]采用了多种传统回归模型，如普通最小二乘线性回归（linear regression，LR）、二次回归和指数回归等，对视野逐点敏感度进行预测，最终发现要达到精确预测未来视野结果，所需的最少视野检查结果数目大约为10。Wen等[21]开发了一种CNN模型5层级联网络（CascadeNet5），该模型仅使用单个视野结果作为输入，即可预测青光眼患眼未来5.5年内的视野逐点敏感度。但是这种模型仅利用单个视野结果中的空间信息，而未充分利用患者过去随访结果中的时间信息。与Wen等[21]的工作不同，Park等[22]的工作更关注于视野随访结果中的时序信息，他们建立了一种递归神经网络（recurrent neural network，RNN）模型，该模型以一系列5个连续的视野结果作为输入，来预测第6个视野结果，该研究预测的是未来视野的总偏差（total deviation，TD）值，而不是逐点敏感度。
对于视野逐点敏感度预测而言，仅基于时间或空间单因素的预测方法可能造成过去视野检查结果中有效信息的浪费，而卷积长短期记忆（convolutional long-short term memory，ConvLSTM）网络是基于长短期记忆（long-short term memory，LSTM）网络的一种改进模型，它可以同时利用数据中的时间信息以及空间信息[23]。已有多项研究证明，在时空序列预测问题上，ConvLSTM模型可表现出相当好的性能[24-26]。因此，为实现更高精度的未来视野预测，本研究采用ConvLSTM网络作为预测模型。本文仅以患者过去视野的随访结果序列作为输入，而无需其它生理结构信息，并将未来0.5～2.0年内视野的逐点敏感度预测作为输出；另外，以样本中包含视野结果的数目作为视野序列长度，对不同长度视野序列对预测结果产生的影响进行了探究，以期找到预测效果最好的视野序列长度。
1 方法
1.1 ConvLSTM网络
本研究旨在联合利用患者过去视野序列的时序特征以及单个视野结果中各点敏感度的空间分布特征，实现对未来视野随访的预测。传统LSTM在输入到状态、状态到状态的转换中采用了全连接结构，无法有效提取空间信息[27]。然而，ConvLSTM将LSTM与卷积操作相结合，将LSTM中的转换结构改进为卷积结构，使得网络具备了对空间结构特征的提取能力，这种方法更适用于兼具长短程时间特性和空间结构特性的视野逐点敏感度预测任务。而且，ConvLSTM在每个时间步都使用相同的卷积核进行卷积操作，实现了参数共享和稀疏连接，相较于全连接结构其参数规模更小，运算效率更高。
ConvLSTM的结构与传统的LSTM基本相似，同样包括遗忘门（F）、输入门（I）、细胞状态（C）和输出门（O）。ConvLSTM的细胞结构示意如图1所示，各部分之间的关系如下：
图1
ConvLSTM细胞结构
Figure1.
ConvLSTM cell structures
（1）遗忘门：它决定了细胞状态中有多少部分被遗忘或保持，其计算过程如式（1）所示：
|
其中，表示卷积；σ(∙)表示S型生长曲线（sigmoid）映射函数；Xt、Ft分别代表第t时间步的输入和遗忘门输出；Ht−1代表第t − 1时间步的隐藏状态输出；Kfx、Kfh分别代表遗忘门在隐藏状态和输入上的卷积核参数，bf为遗忘门对应的偏置参数。
（2）输入门：它决定了当前时间步的新信息有多少部分用于更新细胞状态，其计算过程如式（2）所示：
|
其中，It代表第t时间步的输入门输出，Kix、Kih分别代表输入门在隐藏状态和输入上的卷积核参数，bi为输入门对应的偏置参数。
（3）细胞状态：它结合遗忘门输出和输入门输出来更新自身，其计算过程如式（3）所示：
|
其中，表示哈达玛（Hadamard）积，tanh(∙)为双曲正切函数；Ct − 1、Ct分别代表第t − 1和第t时间步的细胞状态；Kcx、Kch分别代表细胞状态在隐藏状态和输入上的卷积核参数，bc为细胞状态对应的偏置参数。
（4）输出门：它决定细胞状态的哪一部分将作为当前时间步的隐藏状态输出，其计算过程如式（4）～式（5）所示：
|
|
其中，Ot、Ht分别代表第t时间步的输出门输出和隐藏状态输出；Kox、Koh分别代表输出门在隐藏状态和输入上的卷积核参数，bo为输出门对应的偏置参数。
神经网络在广泛的具有挑战性的预测问题上取得成功的原因通常归结于其深度[28]。因此，本研究采用堆叠式ConvLSTM模块，以增加网络深度，从而获得视野序列在不同时间尺度上的特征表示。本研究所采用模型的网络结构包括1个堆叠式ConvLSTM模块和1个卷积层。其中，堆叠式ConvLSTM模块仅包含2个ConvLSTM单元，并且为避免梯度消失，在2个ConvLSTM单元之间添加了1个批标准化（batch normalization，BN）层，对前1个ConvLSTM单元所提取特征的分布进行调整。该模型结构图如图2所示，其中，nVFs表示长度为n的视野序列输入，VFn+1表示模型的预测输出为第n+1个视野结果。同时，n也与ConvLSTM的时间步数对应。由于以时序图表示的单个视野结果尺寸较小，选择较小的卷积核尺寸更适合用于提取视野特征，故设置2个ConvLSTM单元的卷积核大小为3×3，数目为32。最后的卷积层亦选择卷积核大小为3×3，数目为1，以减少通道数目，进而输出视野预测结果。
图2
ConvLSTM网络模型
Figure2.
ConvLSTM network model
1.2 对比模型
为了验证ConvLSTM模型在视野逐点敏感度预测方面的优势，本研究首先构建了4个基线模型，以提供基线预测性能。
本研究选取了2种高斯分布模型（Gaussian distribution model，GDM）和2种LR方法作为基线模型用于后续模型性能比对。其中，第一种GDM模型（记为GDM1）是基于Heijl等[29]工作中的视野进展速率的平均值和标准差。具体而言，本研究使用早期明显青光眼的视野敏感度进展速率数据，即均值为−0.36 dB/a，标准差为0.60 dB/a，将视野序列中逐视野预测结果取平均作为未来视野的最终预测。第二种GDM模型（记为GDM2）则是基于由训练集中统计得到的视野进展速率的平均值和标准差，之后按与GDM1相同步骤进行视野预测[21]。第一种LR模型（记为LR1）针对于患者群体，对视野各点敏感度进行LR操作，得到一组回归系数，以之预测视野。第二种LR模型（记为LR2）则针对于患者个体的逐点敏感度，以患者之前的随访数据预测其未来视野[22]。然后，本研究以ConvLSTM模型结构为基础，构建CNN模型和LSTM模型作为对比模型。CNN模型以单个视野结果作为输入，仅利用视野逐点敏感度的空间位置信息对未来视野进行预测。而LSTM模型则以视野序列为输入，利用序列中的时序信息预测未来视野。CNN模型的结构与ConvLSTM模型基本一致，它由3个卷积层组成，前2个卷积层的卷积核大小、卷积核数目等参数设置与ConvLSTM单元的参数设置基本相同，最后一层卷积层与ConvLSTM模型的最后一层卷积层完全相同。而LSTM模型中同样仅包含2个LSTM单元，并且它以全连接层代替最后的卷积层。
2 实验
本研究使用的开源数据集为来自华盛顿大学汉弗瑞视野分析仪的视野测试数据集（a dataset of perimetry tests from the Humphrey field analyzer at the University of Washington，UWHVF）[30]，在此基础上开展研究工作，以ConvLSTM模型与其他几种对比方法进行预测性能验证与比较。
2.1 数据集
UWHVF数据集中所有患者都接受了视野测量与评估软件汉弗瑞视野分析仪24-2程序模式（Humphrey Field Analyzer，24-2 Program）（Carl Zeiss Meditec Inc.，美国）的视野分析，使用瑞典交互式阈值算法（标准或快速）或全阈值策略执行。数据集包含了3 871位患者的7 428只眼睛的视野数据，共计28 943个视野，所有视野数据均采集自1998年—2018年。每只眼睛的每个视野在时间上进行了对齐，第1个视野被视为每只眼睛的基线视野。对所有纳入数据集的患者进行统计，以“中值[四分位范围]”形式显示：年龄为64[54, 73] 岁；随访的时间为2.49[1.11, 5.03] 年；基线视野的平均总偏差（mean total deviation，MTD）为−4.51[−8.01, −2.65] dB；模式标准差（pattern standard deviation，PSD）为2.41[1.70, 5.34] dB。UWHVF数据集中患者视野示例如图3所示，逐点敏感度精确到百分位，灰度浅深与敏感度大小对应[30]。目前在UWHVF数据集上开展的视野相关的研究较少，尤其在视野敏感度预测方面仍无相关工作发表。
图3
青光眼患者视野与正常视野对比示例
Figure3.
An example of patient's visual field with glaucoma compared to normal person's visual field
2.2 数据预处理
首先，在患者水平上，本研究将UWHVF数据集划分为两部分：80%作为训练集，用于训练和验证阶段；20%作为测试集。两个数据集中患者的统计信息如表1所示（表格仅统计信息完整的数据），可以发现，除性别信息部分缺失导致的差异以外，两个数据集中的数据分布基本相近。为了探究不同长度的视野序列对未来视野的预测性能，本研究将两个集合当中的视野按照时间间距与序列长度分别做了筛选，以得到不同长度的视野序列数据集。具体而言，筛选的标准为：在同一视野序列样本中，相邻视野检查的时间间隔必须在0.5～2.0年内。此区间的选定是考虑到青光眼视野进展较缓慢，时间间距小于0.5年视野变化可能不明显，不利于捕捉其变化规律，而当时间间距大于2.0年，则可认为是视野检测失访，不宜纳入时序研究。此外，为了尽可能多地获得样本，允许视野序列部分重复。例如，某患者的一只眼睛有7个视野测试结果，分别在33.75、34.21、35.00、36.52、37.43、38.65、41.02岁，现可将其分为2个序列长度为4的视野样本：（34.21、35.00、36.52、37.43岁）和（35.00、36.52、37.43、38.65岁），或者1个序列长度为5的样本：（34.21、35.00、36.52、37.43、38.65岁）。并且，部分样本的视野序列中包含全零视野，说明该患者已经造成不可逆致盲，无需纳入研究而予以剔除。
表1
训练集与测试集统计信息
Table1.
Statistics information of training and testing sets
最终，训练集和测试集中样本数目统计如表2所示。由于当序列长度大于7时，数据集中符合筛选条件的样本数目过少，故研究仅针对长度为3～7的视野序列样本（3VFs～7VFs）。对于序列长度为n+1的数据集，实验将以前n个视野结果作为输入序列预测第n+1个视野，记为nVFs→VFn+1，其中n=2, 3 …, 6。
表2
不同长度视野序列数据数目统计
Table2.
Statistics on the number of visual field series data of different lengths
在视野测量与评估软件汉弗瑞视野分析仪24-2程序模式（Carl Zeiss Meditec Inc.，美国）测量下，共获得54个视野测试点。剔除2个生理盲点后，仅使用剩下的52个测试点。每个视野被表示为8×9的矩阵，以保持其空间关系。对所有输入特征采用零均值化方法，将其标准化成均值为0，标准差为1。
2.3 实验参数与评价指标
研究实验平台采用64位Windows10操作系统（Microsoft Corp.，美国），编程语言为Python 3.6.0（Python Software Foundation，美国），而Python中的深度学习算法在TensorFlow 1.7.0（Google Inc.，美国）上使用开源神经网络平台Keras 2.2.0（Francois Chollet et al.，美国）实现。训练深度学习模型的迭代次数均设置为25，批大小为32，选择自适应矩估计（adaptive moment estimation，Adam）作为优化器，损失函数均采用均方误差（mean squared error，MSE）损失。采用5折交叉验证方式对这些模型进行训练，选用交叉验证中表现最优的模型参与后续在测试集上的预测性能对比。
预测模型的评价标准采用平均绝对误差（mean absolute error，MAE）、均方根误差（root mean squared error，RMSE）和决定系数R2，以避免不同评价标准本身存在的局限性。它们对应的计算如式(6)～式(8)所示：
|
|
|
其中，为第i个视野点的真实敏感度，为模型预测第i个视野点的敏感度，为真实敏感度均值，m为预测逐点敏感度的数目。
3 实验结果及分析
3.1 模型对比结果
为了分析ConvLSTM模型在视野预测任务中的表现，使用相同的数据集训练6种对比模型以及来自文献[21]和文献[22]的2种算法模型。其中，文献[22]提出的算法模型是一种RNN模型，因此无需修改模型结构即可接受不同长度视野序列输入。随后，利用各模型在测试阶段的预测结果计算不同的评价指标。结果统计如表3所示，最优值以加粗形式表示。
表3
不同模型预测性能对比
Table3.
Comparison of the prediction performance of different models
通过表3中各个模型的对比结果可以发现：① 从3项评价指标来看，5种深度学习方法的性能显著优于各基线模型，表明它们对视野敏感度变化规律的刻画能力更强，其中，性能最好的是ConvLSTM模型。与仅考虑视野空间或时间单因素信息的CNN和LSTM模型相比，ConvLSTM模型通过联合考虑两种信息，有效提高了视野预测的精度；② 两种GDM模型的性能基本优于两种LR模型，说明使用线性方式难以准确刻画视野敏感度的变化，因此在表现上不如基于数据分布的方法；③ GDM1模型相较于GDM2模型表现更佳，这表明仅依靠训练集得到的统计特征不足以全面描述患者群体特征，因此GDM1在测试集上的预测效果更好；④ LR2模型的表现优于LR1模型，后者的决定系数皆小于0.1，这表明不同患者个体的视野变化差异较大，基于患者整体的LR方法并不适用于视野预测。⑤ 来自文献[21]的CascadeNet5模型与构建的CNN模型相比性能更好，这可能是因为前者的模型结构更加复杂，并且在输入中加入了年龄信息。然而，ConvLSTM模型参数规模远小于CascadeNet5模型，却取得了更好的预测结果，说明时空特征的考虑比增加模型复杂度更加重要。⑥ 来自文献[22]的RNN模型在5种深度学习模型中表现最差，这可能是由于它最初设计用于预测视野的总偏差值而非逐点敏感度。此外，在输入的视野序列数据中人为加入了空白数据，增加了噪声，进一步提升了模型的学习难度。
对于不同缺损程度的视野序列，各深度学习模型的预测性能可能会受到不同程度的影响。本研究将5种深度学习模型应用于不同缺损程度的视野序列进行预测，如图4所示为代表性预测示例，其中包含预测输出值及其误差绝对值两种图像，缺损程度由平均总偏差确定。从图4中示例可以看出，ConvLSTM模型受到视野缺损程度影响相对其他模型更小，主要的预测偏差出现在敏感度数值“分界线”上。进一步地，对ConvLSTM模型在各测试集上的预测结果按照平均总偏差值分段进行统计，采用RMSE指标，并以箱线图方式表达，如图5所示。从图5中的对比结果可以看出，模型对于平均总偏差值较大的测试样本，预测结果基本较好且相对稳定。这表明，对于视野缺损严重的患者，模型能更好地预测其视野敏感度变化，有助于医生对病情发展进行及时干预。此外，由于各测试集中包含的样本数目差距较大，箱线图中“异常值”数目也与之相关，3VFs测试集上的“异常值”最多，7VFs测试集上最少。7VFs测试集中平均总偏差值为“>18 dB”的样本RMSE分布较分散，也是由于样本数目较少所导致。
图4
各模型对不同缺损程度视野的预测代表性示例
Figure4.
Representative examples of predicted visual fields with different defect degrees by different models
图5
ConvLSTM模型在不同序列长度测试集上的预测结果统计（按平均总偏差分段）
Figure5.
Statistics of prediction results of ConvLSTM model on different series length testing sets (segmented by mean total deviation)
3.2 不同序列长度结果
对于时序深度学习模型而言，其输入时间序列的时间步数可能影响最终的预测性能。为了探究不同长度的视野序列输入的影响，本研究对ConvLSTM模型在各测试集上的预测性能进行了比较分析。
由表3中的统计结果可以看出，视野序列长度与模型预测性能之间存在关系。对于ConvLSTM模型而言，结果显示视野序列长度并不是越长或越短预测性能越好，而是长度为3（即3VFs→VF4）的视野序列预测效果最佳。将LSTM和ConvLSTM两种时序深度学习模型的评价指标差值以折线图表示，如图6所示，其中Δ表示差值。可以看到，随着输入视野序列长度变化，MAE、RMSE和R2等指标的差值不断增大，说明二者性能差距在逐步拉开，而且LSTM模型预测指标原值是一直下降的，这表明视野序列中的时序信息的作用在减弱，而空间信息的作用则有所增强。原因在于，随着视野序列长度的增加，视野检测时间跨度增大，带来的时间冗余信息随之不断增多，因此空间信息在视野预测中所占比重上升。在6VFs→VF7的预测中，CNN模型和文献[21]模型预测结果的MAE值反而小于ConvLSTM模型，而其余指标均与后者相差不大，这说明在时间冗余信息的影响下，空间信息所起的作用在逐步凸显。
图6
LSTM和ConvLSTM模型的预测性能差距
Figure6.
Predictive performance gap between LSTM and ConvLSTM models
4 结语
目前，针对视野逐点敏感度预测方面开展的工作还很有限，而现有研究方法仍存在诸多不足，如预测精度低、视野序列中的时空信息利用不充分等。本研究采用ConvLSTM模型对未来0.5～2.0年内视野的逐点敏感度进行预测，以青光眼患者既往视野随访结果作为输入序列，而无需利用其他眼科特检检查或眼部解剖特征，实现了较高精度的预测。主要研究结论如下：① 相较于作为基线的传统模型，深度学习模型在刻画视野逐点敏感度变化方面的能力更强；② ConvLSTM模型充分考虑了视野敏感度随时间变化规律以及各位置点之间的空间关系，因此其预测精度高于CNN模型和LSTM模型；③ 输入预测模型的视野序列长度并非越长越好，时间序列过长可能带来更多冗余信息，实验发现适宜于ConvLSTM模型的输入序列长度为3，即采用过去1.5～6.0年内的3次视野检查结果预测效果更好。
当然，本研究在数据集选用方面存在一定局限性，由于目前开源的视野数据少，研究仅采用了UWHVF数据集，而缺少了其他人群的视野数据。同时，本研究工作是统计学意义上的改进，对于其能否在临床上进行应用，需要结合临床进一步评估，这也将是后续研究需要开展的工作。
重要声明
利益冲突声明：本文全体作者均声明不存在利益冲突。
作者贡献声明：王握负责模型算法的设计和论文的撰写；郑秀娟负责指导方法理论和算法设计，以及指导论文的撰写；吕智清负责资料收集和论文修改；李妮和陈俊负责指导医学基础和分析方法。
