一种基于跨尺度频谱路由的CBCT牙齿分割与牙位编号方法及系统

本发明属于医学图像处理领域，具体涉及一种基于跨尺度频谱路由的cbct牙齿分割与牙位编号方法及系统。

背景技术：

1、cbct能够以较低辐射剂量提供口腔颌面结构的高分辨率三维图像。对每颗牙齿进行准确分割并赋予正确的牙位编号，是正畸方案设计、种植体规划、手术模拟和数字化牙科cad系统的重要基础。该任务需要同时识别牙齿实体、描绘牙冠及牙根边界，并区分牙弓左右侧和上下颌对应牙位。

2、早期方案通常将二维检测框架扩展至三维牙科数据。例如，toothnet采用mask r-cnn检测牙齿边缘；粗到细的v-net方案通过逐级分辨率细化边界；多任务全卷积网络同时预测牙齿区域和表面概率图，再通过标记控制分水岭分离相邻牙齿；基于质心和高斯解耦的方法利用几何中心分开接触牙齿；层级形态学方案还引入候选区域和形状先验。

3、与本技术较接近的临床系统多采用级联流程：首先检测牙齿或牙弓区域，其次进行语义或实例分割，最后通过独立标签模块或启发式规则分配牙位编号。toothseg组合语义分割头与实例分割头，并通过自校正机制处理牙齿合并、分裂及标签重排；oralseg在unetr骨干中结合swin transformer和空间mamba模块，以获得体数据的局部与全局依赖。对于原本只执行语义分割的网络，还需要通过连通域分析拆分接触实例，并在预测实例内采用多数投票分配牙位编号标签。

4、在编码器方面，u-mamba和segmamba等方案用状态空间模型替代或补充自注意力，以线性复杂度建模长程依赖；在跳跃连接方面，unet++采用稠密嵌套连接，attention u-net采用门控抑制无关特征，混合专家方案则把多个编码尺度重采样后交由不同专家处理；在频域方面，傅里叶和小波分解能够提供多分辨率谱信息，但尚缺少将多类频谱描述符与跨尺度跳跃连接路由联合优化的统一牙齿分割方案。

5、现有技术存在如下缺点：

6、1. 多阶段误差传播。区域检测、分割和牙位编号由独立模块依次完成时，前一级漏检、粘连或边界偏差会传递至后一级，且各模块需要分别训练、调参和部署。

7、2. 细粒度边界易丢失。相邻牙齿间隙可能仅有一个体素，牙根与牙槽骨灰度又高度相近；固定的一对一跳跃连接难以根据不同解码层需求同时利用浅层边界细节和深层语义。

8、3. 长程空间关系不足。缺牙、少牙、修复体和个体解剖差异会削弱局部外观线索，牙位编号所需的牙弓全局位置关系难以仅靠局部卷积可靠推断。

9、4. 跨尺度融合缺乏自适应性。简单拼接、相加或注意力门控通常对所有样本采用相似的尺度组合，多个尺度反复重采样还会带来信息冗余和计算负担。

10、5. 稀疏专家可能退化。采用top-k专家路由时，门控容易长期集中于少数专家，导致其余专家训练不足，难以形成互补的空间特征处理能力。

11、因此，本技术提案要解决的技术问题是：在单一端到端网络内完成三维牙齿实例分割与牙位编号；使不同解码层能够按输入自适应获取跨尺度信息；利用全局频率上下文辅助专家选择；保持相邻牙齿和牙根边界；并通过稳定初始化和负载均衡提高训练可靠性。

技术实现思路

1、本发明的目的在于解决上述背景技术问题，提供一种基于跨尺度频谱路由的cbct牙齿分割与牙位编号方法及系统。

2、为实现上述目的，本发明的技术方案是：一种基于跨尺度频谱路由的cbct牙齿分割与牙位编号方法，包括：

3、获取锥形束计算机断层扫描cbct体数据，并为所述cbct体数据附加三个绝对坐标通道，形成网络输入；

4、将所述网络输入送入一u形编码器-解码器网络，所述u形编码器-解码器网络包括混合cnn-mamba编码器、跨尺度频谱路由器cssr、频率上下文聚合器fca、五级解码器、五个深监督分割头以及一个辅助边界头；

5、通过所述混合cnn-mamba编码器提取六级特征金字塔{e0，…，e5}，其中自第2级起插入三方向mamba编码块；

6、通过所述cssr对{e0，…，e4}进行跨尺度自适应重组，得到调制后的跳跃特征；

7、通过所述fca从瓶颈特征e5及深层融合特征提取傅里叶描述符、db2小波描述符和多尺度墨西哥帽描述符，并以film方式调制深层融合特征；

8、通过所述五级解码器进行逐级上采样重建，每一级将上采样特征与对应的cssr调制跳跃特征拼接后送入卷积块，得到五级解码输出；

9、通过所述五个深监督分割头对五级解码输出分别进行33类体素预测，所述33类包括背景类及32颗恒牙类别，从而直接输出逐牙齿实例及其牙位编号；

10、通过所述辅助边界头对最精细解码输出进行单通道边界图预测。

11、进一步的，所述三方向mamba编码块对输入的处理包括以下子步骤：

12、对x进行组归一化和深度可分离3×3×3卷积，得到；

13、分别沿深度d、高度h和宽度w三个方向将重排为序列，经层归一化和mamba状态空间层处理后恢复空间形状，得到三个方向的状态空间输出；

14、将三个方向的状态空间输出按通道平均，经1×1×1卷积投影后与输入x残差相加，得到中间特征xmid：

15、

16、将xmid送入由两层1×1卷积和gelu构成的前馈网络，再次与xmid残差相加，得到编码块输出：

17、。

18、进一步的，所述cssr包括组a和组b两个并行工作的混合专家组，其中：

19、组a面向{e0，e1，e2}并以e1为锚点，各特征先投影，再通过组归一化、relu和带可学习残差权重的深度卷积进行空间细化，各特征重采样至锚点分辨率后拼接融合，进入4个空间专家；组a设置3个独立门控网络，分别产生返回e0、e1和e2的输出；

20、组b面向{e2，e3，e4}并以e3为锚点，采用比组a更大的投影维度，使用原始未调制的e2，并设置2个独立门控网络；

21、组a对e2产生的残差与组b并行计算，不串联传播；

22、cssr对每一跳跃特征ei的调制公式为：

23、

24、其中oi为门控专家输出，proji为映射回原通道数的两层瓶颈结构，αi从0开始学习。

25、进一步的，每个空间专家的内部结构包括：输入1×1卷积、组归一化、relu、3×3×3深度卷积、输出1×1卷积、组归一化、relu；其中输入投影将通道数加倍，输出投影恢复通道数；

26、每个门控网络对全局平均值、全局最大值和空间标准差进行拼接，经小型多层感知器及softmax获得专家权重：

27、

28、其中，g(z)表示门控网络针对输入特征z生成的专家权重向量，其长度等于专家数n；z表示送入门控网络的融合特征；mlp表示多层感知器；gap(z)表示对z在空间维度执行全局平均池化；gmp(z)表示对z在空间维度执行全局最大池化；stdspatial(z)表示z在空间维度的标准差；方括号内的分号表示特征拼接；softmax表示将门控得分归一化为专家权重。

29、每个门控仅激活权重最高的k=2个专家，并对这两个权重再次softmax归一化；组a的3个门控与组b的2个门控合计m=5个独立门控网络，每个门控对应n=4个专家;

30、按训练批次累计各门控对各专家的选择次数um，并采用使用次数的平方变异系数作为负载均衡项：

31、

32、当各专家使用量近似均衡时上述损失接近0；单一专家占主导时上述损失增大；对于选择次数为0的专家，周期性抽取随机批次元素进行诊断前向传播，并对其输出施加l2惩罚。

33、进一步的，所述fca包括三个并行频域分支：

34、傅里叶分支：从瓶颈e5降维后执行三维实数快速傅里叶变换，令a = log(1 + |rfft(e5)|)，从a提取全频网格通道均值、最大值、低频区均值和高频区均值，低频区取各轴频率网格前四分之一，高频区取其补集，拼接后得到傅里叶描述符vfft；

35、rfft(·)表示三维实数快速傅里叶变换；|·|表示复频谱幅值；

36、db2小波分支：采用固定daubechies-2三维外积滤波器，以冻结权重的分组卷积将组b融合特征分解为8个子带，对每个子带计算平均绝对值和均方根后拼接；

37、墨西哥帽分支：采用三个尺度递增的改进ricker小波核

38、，对组b融合特征执行固定权重的深度三维卷积，每个核做零均值和l2归一化，按尺度提取平均绝对值、均方根和最大绝对值；

39、三路描述符分别经layernorm-linear-relu-linear-relu投影到公共维度后拼接，再经layernorm-linear-relu-linear产生逐通道参数和，以film方式调制组b融合特征fb：

40、

41、fca融合网络最后一层采用零初始化，使训练开始时==0，fca最初为恒等调制。

42、进一步的，所述五级解码器的解码过程为：

43、以d5=e5作为起始解码特征；

44、第l级解码输出为：

45、

46、五个解码输出均连接一个深监督分割头，产生33类体素预测；

47、最精细解码输出d0还连接轻量级边界头，预测单通道边界图；

48、真实边界b由相邻体素类别是否不同确定，且两相邻体素至少一个属于前景牙齿，不对边界执行膨胀。

49、进一步的，总损失函数为：

50、

51、其中β(t)采用先预热后余弦衰减调度；

52、分割损失对五个深监督输出同时计算，由困难体素topk交叉熵与tversky损失等权组合：

53、

54、表示五个深监督输出上的分割损失；表示仅对损失最高的预定比例体素计算的困难体素交叉熵；表示仅对当前批次实际出现的牙齿类别计算的tversky损失；和分别表示上述两项损失的权重，二者等权；

55、困难体素topk交叉熵仅计算损失最高的预定比例体素；tversky损失只统计当前批次实际出现的牙齿类别：

56、

57、其中β＞α，以提高对假阴性的惩罚；

58、边界损失为边界预测与真实边界之间的二元交叉熵和dice损失之和：

59、

60、lboundary表示边界监督损失；b̂表示辅助边界头输出的预测边界概率图；b表示由真实分割标签生成的真实边界图；lbce(·,·)表示二元交叉熵损失；ldice(·,·)表示dice损失。

61、进一步的，训练阶段采用以下配置：

62、所有cbct扫描重采样至0.3mm各向同性间距，训练块尺寸为128×160×160；

63、使用4块32gb gpu、总批量12，训练1000个epoch；

64、优化器为带nesterov加速的sgd，动量0.99，权重衰减3×10-5，初始学习率0.01并采用多项式衰减；

65、启用fp16自动混合精度和torch.compile；

66、推理时采用0.5重叠滑窗与高斯融合，不使用测试时增强。

67、进一步的，所述网络输入在送入混合cnn-mamba编码器之前，以nnu-net的自配置阶段设计为基础进行预处理；所述混合cnn-mamba编码器的六级卷积均采用3×3×3卷积核，卷积后采用实例归一化和leakyrelu。

68、本发明还提供了一种基于跨尺度频谱路由的cbct牙齿分割与牙位编号系统，执行如上述所述方法，包括：

69、图像获取模块，被配置为获取cbct体数据，并为所述cbct体数据附加三个绝对坐标通道，形成网络输入；

70、混合编码模块，被配置为通过混合cnn-mamba编码器提取六级特征金字塔{e0，…，e5}，其中自第2级起插入三方向mamba编码块；

71、跨尺度频谱路由模块，被配置为通过cssr对{e0，…，e4}进行跨尺度自适应重组，得到调制后的跳跃特征；

72、频率上下文聚合模块，被配置为通过fca从瓶颈特征e5及深层融合特征提取傅里叶描述符、db2小波描述符和多尺度墨西哥帽描述符，并以film方式调制深层融合特征；

73、解码与预测模块，被配置为通过五级解码器进行逐级上采样重建，每一级将上采样特征与对应的cssr调制跳跃特征拼接后送入卷积块，得到五级解码输出；通过五个深监督分割头对五级解码输出分别进行33类体素预测，并通过辅助边界头对最精细解码输出进行单通道边界图预测；

74、其中，所述33类包括背景类及32颗恒牙类别，网络输出直接给出逐牙齿实例及其牙位编号。

75、相较于现有技术，本发明具有以下有益效果：

76、1. 减少级联误差。由于分割和牙位编号类别在同一网络中联合输出，检测漏失或后处理规则不会作为独立阶段继续放大误差，部署链路也更简洁。

77、2. 提高边界与语义的协同利用。双组cssr使每个解码尺度根据样本从多个编码尺度选择信息，浅层细节与深层语义不再受固定跳跃连接限制。

78、3. 增强全局位置判别。fca把傅里叶全局谱、小波子带能量和多尺度局部频率响应注入门控，为缺少局部解剖参照时的牙位编号判断提供额外上下文。

79、4. 保持训练稳定。零初始化使cssr和fca在训练初始不破坏骨干网络表征，负载均衡损失减少专家长期闲置。

80、5. 改善困难体素和牙齿边界。topk交叉熵集中优化难分体素，tversky提高对漏分牙齿的惩罚，边界头直接约束相邻类别交界。

81、6. 实验效果明确。在toothfairy 3测试集上，本技术方法实例f1为0.9783、全景dice为0.9280；要求牙位编号标签正确的对象f1为0.9650、全景质量dice为0.9160，均高于所列六种对比方法。