2001年5月30日，随着中国东方航空从美国芝加哥至中国上海的MU588号航班顺利着陆，标志着中国民航“零的突破”——首次实现极地飞行。此次极地航路试飞成功，使中国与美国东海岸之间1.4万km航程缩短了1 300多km^[1]。所谓极地航线，是指穿越极地区域的飞行。极地航线由于采用大圆航线策略，可显著缩短航程，在满足航路导航、通信和气象条件的前提下，有助于减少燃油消耗和碳排放^[2]。然而，极地航线在带来经济效益的同时，也引发了对航空机组人员职业健康风险的高度关注^[3]。由于地球磁场在极地区域相对较弱，宇宙高能粒子更易穿透大气层进入对流层，造成大气中次级粒子的增强，从而使极地上空的宇宙辐射水平显著高于中低纬度地区^[4]。在高空飞行过程中，尤其是巡航高度通常在10~12 km的民航客机^[5]，航空机组人员所受辐射主要来自宇宙射线，这种辐射强度受飞行高度、纬度、飞行时间和太阳活动等因素的共同影响^[6]。2023年11月，韩国一名航空机组人员罹患胃癌死亡，被韩国劳动福祉工团和职业病(业务疾病)裁定委员会认为与其工作存在一定程度的因果关系，认定为工伤。这一事件被国内媒体广泛报道，引起了社会对航空机组人员宇宙射线照射及其健康危害与控制的关注^[7]。有研究表明，机组人员职业性辐射暴露可能与白血病、乳腺癌、皮肤癌等疾病的发生风险升高有关^[8-10]。此外，宇宙辐射还可能导致生殖细胞损伤、胚胎发育异常等问题，尤其对处于妊娠期的飞行人员构成较高风险^[11]。因此，国际放射防护委员会(ICRP)强调，应加强对航空辐射暴露的评估与监测，合理控制机组人员的年累积有效剂量^{[4, 12]}。

1.航线采集：通过Flightaware网站查询航路穿越北纬78°以上区域的民航航班，记录其始发地与目的地，剔除往返均使用极地航路的航班。此外，针对如有多家航空公司运营同一条航线，或同一航空公司同一航线每日多班次的情况，由于航路相同，按同一条航线处理。经分析，共筛选出26对52条航线，覆盖北美—东亚、东亚—欧洲等典型跨极区航段，其往返航程中均包含一程穿越北纬78°以上区域的极地航线及一程不经极地的对照航线。考虑商业数据使用限制，本研究不逐一列出每条航线的实时航路点。查询并记录上述52条航线的飞行高度、飞行时间、飞行实时经纬度等飞行轨迹数据，完成航线数据采集。

2.航空机组人员所受有效剂量计算：根据中国科学院国家空间科学中心太阳活动周监测数据(http://www.sepc.ac.cn/dailyForecast_chn.php)，第24个太阳活动周期为2008—2019年。在本太阳活动周期内太阳活动水平具有典型的“强-中-弱”三阶段特征，其中2014年处于太阳活动极大期、粒子调制最强；2017年位于极大期向极小期过渡的中间阶段，可反映中等调制背景；2019年则处于太阳活动极小期、宇宙射线通量最高。本研究选取2014年4月、2017年2月和2019年12月作为航线剂量计算的代表时间点。选取这3个时间点既能覆盖整个太阳活动周期的典型状态，又可在保证计算可行性的前提下反映太阳调制对航线宇宙辐射剂量的年度差异，从而避免因单一年份计算导致的外推偏差。

采用CARI-7A计算不同太阳活动水平下所选择航线沿着指定路径飞行的航空机组人员所受有效剂量。CARI-7A为美国联邦航空管理局(Federal Aviation Administration，FAA)提供的航空宇宙辐射剂量评估工具，其计算架构具有多项可配置参数以适配不同研究需求。程序在初级宇宙射线谱的选择上支持多种模型，例如内置的国际标准化组织(International Organization for Standardization，ISO)银河宇宙射线模型及其他经验/理论谱，在大气输运计算上既可采用叠加近似以简化重离子级联为多个独立核子的问题，也可在需要时采用更细致的处理；在剂量换算上可选用不同版本的辐射权重因子(包括但不限于ICRP 60报告^[12]、ICRP 103报告^[13]等)；在轨迹输入方面既支持以出发/到达机场代码为输入的默认大圆航线计算，也允许用户按经纬度自行定义航路点(包含高度)以精确描述实际飞行路径。上述可选项使CARI-7A在保持计算效率的同时，可根据研究对象与对照要求在精度与可比性之间进行权衡。

基于上述灵活性与本研究目标，本研究在CARI-7A中作出如下明确配置以保证物理一致性与可比性：①选用ICRP 103号报告^[13]中推荐的辐射权重因子。②选用ISO推荐的银河宇宙射线(galactic cosmic ray，GCR)初级谱模型^[14]，一方面ISO模型为多项国际标准和既有航空剂量研究广泛采用的基线谱，便于与国内外前人结果进行比较；另一方面，本研究跨越不同太阳活动阶段(2014、2017、2019年)且关注长期背景剂量差异，ISO模型在对GCR基线与随太阳调制变化的表征上提供了稳定且保守的基准，有利于减小因谱模型选择带来的系统偏差。③在大气级联处理上采用叠加近似法，该方法在航空剂量计算领域被广泛采用，能够在保留对重离子贡献物理意义的同时显著降低计算复杂度，从而在对大量航班进行批量计算时兼顾效率与物理可解释性，同时，叠加近似对于常规巡航高度(10~12 km)及本研究所涉及的粒子能谱区间具有良好近似精度，因此适合本研究的对比目的。④计算输入采用自行定义的航路点(经纬度与巡航高度)而非仅以机场代码为输入的默认大圆航线。⑤输出的剂量指标选择有效剂量。有效剂量是按照ICRP推荐的辐射防护量定义的综合性剂量学指标，其数值由多种粒子(包括中子、γ射线、质子、电子/正电子、μ介子及少量重离子等)结合各自的辐射权重因子共同换算得到，能够整体反映飞行人员在复杂宇宙射线场中的受照风险。CARI-7A在大气级联输运计算中会模拟不同类型的初级和次级粒子通量，并在计算有效剂量时将所有粒子贡献自动整合为一个总剂量值^[15]，因此本研究使用的有效剂量结果已包含主要辐射组分(特别是中子与光子)的剂量贡献。基于以上配置，本研究在保证与既往文献可比性的同时，最大化地保留了对极地与非极地路径物理差异的解析能力，从而为后续线性混合效应模型中关于“是否为极地航线”主效应的量化提供了稳健且可重复的剂量输入数据。

本研究中所有剂量结果均指单次航班的航空机组人员所受有效剂量，即沿该航班实际飞行路径从起飞至落地的航空机组人员所受有效剂量，不涉及年度累积剂量或多次飞行的剂量合并。

3.数据分析：考虑到太阳活动水平、飞行时间及巡航高度等因素对航空辐射剂量的影响，本研究引入太阳射电流量指数(F_10.7)、地磁活动指数(K_p)、飞行时间(T)与平均巡航高度(H_mean)作为协变量纳入分析。F_10.7指数是指太阳在10.7 cm波段(相当于2 800 MHz频率)发出的射电流量，其数值以每平方米每赫兹的瓦数(sfu，1 sfu = 10^-22 W·m^-2·Hz^-1)表示，是衡量太阳活动强度的重要物理指标；F_10.7指数越高，表明太阳活动越强，对宇宙射线的调制作用增强，从而影响大气层中高能粒子的通量；Kp指数则反映地磁活动强度，范围通常在0至9之间(无量纲)，数值越高代表地磁扰动越强，可能降低地磁场对宇宙射线的屏蔽能力，进而影响辐射剂量。为减少每日指数波动对剂量估算结果的影响，本研究采用F_10.7和Kp指数的月平均值进行统一计算。通过美国国家海洋和大气管理局(National Oceanic and Atmospheric Administration，NOAA)网站获得F_10.7和Kp指数历史数据。

本研究采用线性混合效应模型(linear mixed-effects model，LMM)评估极地航线对航空机组人员所受有效剂量的影响，选择该模型的主要理由在于研究数据具有典型的分层结构：每对往返航班包含一条极地航线与一条非极地航线，属于同一航线对的两条记录存在潜在相关性，而不同航线对之间相互独立。若采用传统线性回归模型，无法处理这种嵌套结构，可能低估标准误并导致假阳性结果。LMM通过引入航班编号作为随机效应，可显式建模同一航线对内部的相关性，保证参数估计的无偏性和方差估计的准确性。分析时将“是否为极地航线”作为主要固定效应，同时纳入H_mean、T、F_10.7及Kp等协变量，以航班编号作为随机效应项，以控制同一航班重复观测带来的相关性。分析前对因变量进行Lilliefors正态性检验和Levene方差齐性检验，并通过方差膨胀因子(VIF)评估多重共线性。模型采用最大似然法(ML)估计固定效应参数，显著性水平设定为α = 0.05，并计算边际R²和条件R²用于衡量固定效应与整体模型的解释力。

1.航线辐射剂量描述性统计：本研究共纳入26对国际航线往返航班数据(共156条航班记录)，覆盖2014年4月、2017年2月和2019年12月3个时间点。每对航线中均包含一程经过北纬78°以上(极地航线)，另一程未经过(非极地航线)，出发地与目的地相反。为便于展示后续模型，表达式中用符号λ代表航线类型(1=极地航线；0=非极地航线)。所有航班均记录了航空机组人员所受有效剂量、T、H_mean、F_10.7和K_p指数。总体结果显示，极地航线航空机组人员所受有效剂量为(61.43±20.18)μSv，高于非极地航线的(47.94±14.90)μSv。极地航线的飞行时间和平均巡航高度略高于非极地航线，但数值差距不大。根据NOAA提供的历史数据，所采用的3个时间点对应的F_10.7月均值和K_p月均值分别为：2014年4月，F_10.7=141.6 sfu，K_p=2.9；2017年2月，F_10.7=79.5 sfu，K_p=3.0；2019年12月，F_10.7=70.2 sfu，K_p=1.7。3个时间点的F_10.7指数显示太阳活动逐年减弱。对应的K_p指数变化幅度相对较小。表 1列出2014年4月、2017年2月及2019年12月3个时间点的航线数据特征，其中每年均涵盖26对往返航线的极地/非极地航空机组人员所受有效剂量的均值及范围、F_10.7与K_p指数，直观反映太阳活动减弱背景下剂量水平的年度趋势。对飞行时间与平均巡航高度的分析显示，极地航线(78条)的飞行时间为(859.58 ± 58.50)min，范围为736~979 min；其平均巡航高度为(10 802 ± 1 078)m，范围为8 458~12 192 m。而非极地航线(78条)的飞行时间为(795.69 ± 69.71)min，范围681~930 min；平均巡航高度为(10 765 ± 802)m，范围8 586~11 895 m。图 1展示了26对航线在不同时间下的航空机组人员所受有效剂量分布情况。结果表明，极地航线对应的航空机组人员所受有效剂量普遍高于非极地航线，反映了飞行纬度对宇宙射线辐射暴露水平的影响。

表 1(Table 1)

表 1 极地与非极地航线航空机组人员所受有效剂量对比 Table 1 Comparison of effective doses received by aircrew on polar and non-polar routes 年份 航线类型 Emean±SD(μSv) 有效剂量范围(μSv) F10.7 (sfu) Kp 2014.04 P 60.50±10.43 41.44~82.37 141.6 2.9 N 44.31±9.69 33.95~63.73 2017.02 P 60.85±10.49 41.67~82.85 79.5 3.0 N 44.41±9.76 34.27~64.02 2019.12 P 63.00±10.86 43.15~85.90 70.2 1.7 N 45.39±10.29 34.42~66.02 注：P.极地航线；N.非极地航线；F10.7. 太阳射电流量指数; Kp.地磁活动指数

表 1 极地与非极地航线航空机组人员所受有效剂量对比 Table 1 Comparison of effective doses received by aircrew on polar and non-polar routes

图 1 不同年份航线航空机组人员所受有效剂量分布图 Figure 1 Distribution of effective doses received by aircrew on different routes in various years

2.线性混合效应模型

(1) 适用性验证：在模型建立前，分别进行了正态性与方差齐性检验。Lilliefors正态性检验结果为P = 0.001，提示因变量略偏离正态，但由于线性混合效应模型对正态性的稳健性较强，仍可继续分析。Levene方差齐性检验P = 0.51，表明数据方差齐性假设成立。此外，对主要自变量进行了多重共线性诊断，结果显示所有变量的方差膨胀因子VIF均＜5(最大值为1.60)，不存在严重共线性问题，满足模型使用前提。

(2) 模型构建与参数估计：本研究构建的线性混合效应模型公式如下：

式中，λ为航线类型(1 =极地航线；0 =非极地航线)；u_i为随机效应，ε_ij为残差项。模型包含156条观测记录，26个随机效应水平。模型拟合采用最大似然估计，相关模型结果列于表 2。线性混合效应模型的随机效应与残差的估计结果如下：模型随机效应的标准差为7.70 μSv，方差为59.33 μSv²，残差的标准差为0.86 μSv，方差为0.73 μSv²。基于方差分量计算的组内相关系数ICC为0.99(99%)，表明组间变异占总变异的99%，数据具有强烈的聚类特征。因此，使用混合效应模型处理数据的非独立性是必要且合适的。模型的边际决定系数Marginal R²和条件决定系数Conditional R²均为0.99，表明模型拟合良好。

表 2(Table 2)

表 2 线性混合效应模型建立结果 Table 2 Results based on established linear mixed effects model 变量 β 标准误 t值 P值 95% CI 截距项(μSv) -83.55 23.60 -3.54 ＜0.05 -130.17，-36.92 航线类型(μSv) 9.61 2.45 3.92 ＜0.05 4.76，14.46 平均飞行高度(μSv·m-1) -0.10 0.15 -0.68 0.50 -0.40，0.19 飞行时间(μSv·min-1) 0.004 0.001 3.17 ＜0.05 0.002，0.007 F10.7指数(μSv·sfu-1) 0.11 0.02 5.94 ＜0.05 0.07，0.15 Kp指数(μSv) -0.02 0.002 -7.99 ＜0.05 -0.02，-0.01

表 2 线性混合效应模型建立结果 Table 2 Results based on established linear mixed effects model

每个回归系数β表示该变量每变化一个单位时对航空机组人员所受有效剂量(μSv)的边际贡献。若为正值，表示变量增大，剂量也增大；反之亦然。模型结果显示，在控制上述协变量的基础上，是否为极地航线对航空机组人员所受有效剂量的影响具有统计学意义(β = 9.61，P＜0.05)，表明极地航线航空机组人员所受有效剂量显著高于非极地航线。此外，飞行时间(β = 0.004，P＜0.05)和F_10.7指数(β = 0.11，P＜0.05)与航空机组人员所受有效剂量亦呈显著正相关，而K_p指数则呈负相关关系(β = -0.02，P＜0.05)。在本研究中，平均巡航高度对航空机组人员所受效剂量的影响差异无统计学意义(P > 0.05)，与既有理论中飞行高度越高，宇宙射线强度越大，剂量越高的认识存在一定差异。造成该结果的原因可能与本研究所使用的数据样本中航班平均巡航高度分布相对集中有关，样本间高度差异较小，难以充分体现高度对辐射剂量的梯度效应。此外，巡航高度的单一性可能弱化了其在多因素分析模型中的解释力。因此，后续研究中可考虑纳入更多不同高度梯度的航班数据，以更准确评估高度因素对航空辐射剂量的影响。同理，本研究将太阳射电流量指数和地磁活动指数作为控制变量纳入模型分析，结果虽提示其对航空机组人员所受有效剂量的影响具有统计学意义，然而由于样本中这两个指标的水平分布较为有限，仅涵盖了3个不同数值，变量变化范围较窄，可能导致模型对其效应的估计缺乏外推性和稳健性，因此相关结果应谨慎解读。此外，F_10.7和K_p作为协变量的引入，主要目的是控制其对宇宙射线水平的干扰作用，以更准确地评估是否为极地航线对航空机组人员所受有效剂量的主效应影响，其本身并非研究关注重点。尽管变量分布存在一定局限性，在当前样本条件下其纳入有助于提升模型的控制能力，整体上不会实质性削弱主效应变量解释力的可信度。

本研究基于26对极地与非极地航线，利用CARI-7A程序估算不同航线路径下航空机组人员职业照射水平，并结合线性混合效应模型对是否为极地航线与剂量水平之间的关系进行了定量分析。结果显示，在控制飞行时间、平均巡航高度、F_10.7和K_p等协变量的前提下，极地航线对辐射剂量具有显著影响，航空机组人员所受有效剂量显著高于非极地航线，进一步印证了纬度对航空宇宙射线暴露水平的关键作用。这一结论与既往研究一致，如拓飞等^[2]通过四元件光致发光(OSL)剂量计和CR-39径迹探测器测定发现极地航线空勤人员的有效剂量高于非极地航线。ICRP亦指出，纬度升高将导致宇宙射线通量显著增强，是航空剂量最重要的地理影响因素之一^[4]。航线路径的选择直接影响纬度与经度路径，事实上已有研究表明航路变更(如战争空域绕飞)会显著改变航空机组人员辐射暴露，甚至在极端情形下增加50%^[16]。因此，基于本研究结果与国际实践，建议采取可行措施以降低或管理极地航线所致的职业暴露风险。

此外，尽管本研究的建模结果强调极地航线本身对航空机组人员所受有效剂量的显著贡献，但模型对高度与时间的分析敏感度可能不足。一方面，由于样本平均飞行高度范围集中，难以探测高度差异对剂量的细微影响；另一方面，太阳高能粒子事件(solar particle events，SPE)可显著增加高纬度航段的辐射暴露，尤其是在发生地面增强事件(ground level enhancement，GLE)时，其影响更为剧烈。但是，本研究并未考虑此类极端情形。例如，有模拟研究显示，在典型GLE期间，纬度60°以上、飞行高度约12 km的区域，有效剂量率可较平常增加近50%^[17]。即便是在未触发GLE的强SPE事件中，常规巡航高度(约12 km) 的有效剂量率亦可能升高至4.5 μSv/h^[18]，显著高于正常背景水平。由此可见，如果研究样本未涵盖这些强太阳事件及其时间空间变化，将可能低估航线实际的剂量极值及其变化幅度。因此，在下一步研究中，将扩大巡航高度与样本跨度，进一步评价巡航高度与飞行时间对航空剂量的调制效应。

另外，本研究旨在量化“极地”与“非极地”航线这一分类下的平均剂量差异，所报告的约13 μSv(即模型系数β=9.61 μSv，与描述性统计差值趋势一致)是基于所选的26对特定航线的计算结果。商业极地飞行只位于北极地区，且受航路审批、空中交通管制等因素影响，实际路径存在差异。因此，这一差异值可作为该类航线辐射风险增量的一个代表性估算，但其具体数值会因实际飞行的精确纬度范围、飞行剖面(高度与时间的动态组合)以及所选取的对照航线而波动。未来研究可通过构建更大规模、更具地理多样性的航线数据库，来进一步验证和细化这一差异的普遍范围。其次，各影响因素(如纬度、太阳活动、地磁活动、飞行高度与时间)对总剂量的贡献并非简单的“有”或“无”，而是存在复杂的交互与份额关系。本研究通过线性混合效应模型，在控制其他协变量的情况下，成功分离并量化了“极地航线”这一综合纬度因素的主效应。然而，要更精细地解析纬度(作为连续变量)、太阳/地磁活动以及飞行参数之间的定量关系与相对重要性，需要采用更复杂的模型(如包含交互项的多变量模型)或基于物理的辐射传输模拟，在更广泛的数据范围内进行系统分析。这将是后续深入研究的重点。

利益冲突  无

作者贡献声明  张景国负责文献调研、调查数据录入、整理、汇总、分析和论文撰写；张雨负责调查和数据采集、复核；周文珊和范胜男参与数据录入与分析；李梦雪参与数据整理与分析；邓君指导研究设计和论文撰写