本發明屬于生物醫學成像，具體涉及一種基于改進磁化模型和先驗信息驅動的磁粒子重建方法。

背景技術：

1、磁納米粒子成像(以下簡稱mpi)是一種高對比度、高靈敏度、高空間分辨率的分子影像技術，它使用超順磁性氧化鐵納米粒子作為示蹤劑，利用粒子的非線性磁化特點顯示疾病早期的分子、蛋白質異常，從而實現了精準診療。mpi的空間分辨率直接決定了成像質量與診療的精準度，因此提出提升空間重建精度的新方法是十分必要的。在mpi中，利用永磁體或通電線圈產生的梯度場來創建無場區(ffr)。然后疊加一個正弦振蕩驅動場，使無磁場區沿著預定軌跡(如利薩如軌跡)掃描視場(fov)。在ffr內，納米粒子可根據驅動磁場調整其磁化，接收線圈捕捉這種非線性響應以產生電壓信號。然后，根據電壓信號來可視化待測區域中的納米濃度分布情況，進而對被測物內部進行直接成像。

2、磁納米粒子成像技術的精準度與系統矩陣的構建方法密切相關，系統矩陣由每個磁性納米粒子在所有可能的空間位置產生的時間信號的傅里葉分量組成，表示測量信號與納米顆粒分布之間的關系。基于sm的mpi重建的準確性和質量在很大程度上取決于sm的獲取方式。目前，構建sm的方法主要有兩種：基于測量的方法和基于模型的方法。在sm測量記錄相應的信號時，必須在fov內的每個體素上重復移動一個增量mnp樣本。為了提高sm測量的質量，通常采用多重平均法，顯著延長了校準時間。此外，當示蹤劑的性質或磁場環境發生變化時，sm需要重新校準，從而帶來大量的人工和時間成本。基于模型的方法通常使用能夠代表成像系統和納米粒子物理特性的數學模型構建sm。這些模型考慮了磁場分布、納米粒子特性和掃描程序等因素，無需進行大量測量。基于朗之萬模型獲得sm的可行性在許多研究中得到了證實，然而，朗之萬模型忽略了spios的弛豫效應和相互作用，這對實際情況中精確模擬spios的磁化強度帶來挑戰。因此，引入合適的模型來準確描述mpi中磁性納米粒子的動態磁化行為對mpi電壓分析和提高mpi重建結果的分辨率至關重要。

3、在理想條件下，基于sm的mpi圖像重建可以被認為是一個線性逆問題，假設濃度分布中粒子間的復雜相互作用最小。各種傳統的正則化解決方案，特別是與迭代方法相關的解決方案，如著名的kaczmarz算法和共軛梯度法(cgs)算法，已經被提出用于mpi重建。此外，還引入了先進的正則化模型，如非負融合lasso技術、包含最小-最大凹面的非凸正則化技術和tv懲罰正則化技術，以增強圖像的清晰度和均勻性。雖然這些方法提高了重建精度，但它們通常需要復雜的矩陣操作和手動參數調整，這可能會影響穩定性和圖像質量。

4、近年來，將基于深度學習的算法應用于mpi重構已成為一種趨勢。深度圖像先驗(dip)方法通過利用在測試期間學習的未經訓練的網絡來增強泛化。另一種方法，即插即用mpi(pp-mpi)方法，涉及預訓練圖像去噪，并將其與split-merge(sm)算法相結合，以實現效率和泛化之間的平衡。然而，這些方法在跨不同任務傳遞優先級時可能面臨性能限制。并且大量標簽數據的來源和模型可解釋性較差也一大障礙。

技術實現思路

1、為了克服上述現有技術的不足，針對mpi重建圖像精度質量不足的問題，本發明提出一種基于改進磁化模型和先驗信息驅動的磁粒子重建方法，處理對象基于磁納米粒子成像技術，以改進的jiles-atherton磁滯函數作為磁納米粒子磁化模型，采用梯度場作為選擇場，正弦變化的均勻磁場作為驅動場，利用接受線圈采集磁化產生的感應電壓，隨后結合磁納米粒子成像過程中的先驗信息(飽和磁矩、粒子濃度、真空磁導率等)，并通過深度展開模型進行磁納米粒子濃度的重建，從而定量進行磁粒子濃度在生物組織體內的可視化展示。該方法可用于求解mpi濃度重建的逆問題，解決了傳統重建算法中迭代過程以及人為調參費時費力的缺點，進行高效且穩定的mpi圖像重建，獲取了高質量的圖像重建結果。

2、為了實現上述目的，本發明采用的技術方案是：

3、一種基于改進磁化模型和先驗信息驅動的磁粒子重建方法，包括以下步驟：

4、步驟一、基于磁納米粒子成像的基礎理論，結合經過改進的jiles-atherton磁化模型，用以獲取系統矩陣(表征磁場與磁化強度關系的矩陣)和粒子分布的測量電壓(反映磁化強度變化導致的感應電壓)；

5、步驟二、考慮磁納米粒子解空間的獨特性質及空間相鄰體素間磁化的內在關聯性，分別設計了兩種先驗信息。一是l1正則化先驗，旨在促進解的稀疏性，有效抑制噪聲干擾，提高重建圖像的質量；二是磁化權重的laplacian正則化先驗，該約束通過描述相鄰粒子空間中的磁化梯度關系，進一步增強了重建結果的局部連貫性和空間一致性；

6、步驟三、基于步驟二中提出的先驗信息，構建了一個雙先驗約束的目標函數。該目標函數旨在平衡數據保真度與解的稀疏性、平滑性之間的關系，為后續的優化求解奠定堅實基礎；

7、步驟四、為了對目標函數進行高效迭代求解，基于非凸優化和迭代逼近理論，構建深度展開模型，以到達對磁納米粒子濃度分布的高保真度重建；

8、步驟五、為了驗證所提方法的實際應用效果，設計合理的損失函數進行模型訓練并用測試數據集進行預測。同時，展示了預測圖像和真實圖像。

9、進一步，步驟一具體方法步驟如下：

10、在磁納米粒子成像(mpi)過程中，系統主要由三部分組成：選擇場、驅動場和接收線圈，首先，選擇場是一種梯度場，產生零場區(ffr)，該區域的粒子由于未達到磁化飽和，從而有產生磁化響應的前提條件，產生零磁場的區域一般分為零磁場點(ffp)和零磁場線(ffl)，本發明基于ffp進行研究。其次，疊加一個時變驅動場使ffr在視場(fov)內移動，該驅動場作為ffr內磁納米粒子的激勵場，使粒子產生磁化響應，接收線圈根據磁通量密度變化而產生感應電壓信號。本發明基于磁學的基本原理結合磁納米粒子的具體特性，采用改進的jiles-atherton磁化模型對磁納米粒子在mpi中的磁化行為進行了精確建模，用以獲取系統矩陣和粒子分布的測量電壓。首先，考慮磁納米粒子與外場之間的復雜相互作用，定義了總有效磁場he為：

11、he＝h+αm??(1)

12、其中，h是外場，α是耦合參數，m是磁納米粒子的磁化強度。其次，應用朗之萬函數作為無滯后磁化函數，將無磁滯磁化man與有效磁場he的關系表示為：

13、

14、其中，是郎之萬參數，ms為磁納米粒子的飽和磁化強度。在mpi中，可逆磁化mrev與疇壁彎曲有關，不可逆磁化mirr與疇壁位移有關。磁納米粒子沿外磁場排列，導致疇壁位移。因此，不可逆磁化強度mirr占主導地位，改進的jiles-atherton磁化模型mja可以寫為：

15、

16、為了簡化模型，可以假設可逆參數c為零。此時，mja＝mirr。mja模型可以寫成描述粒子在外加磁場h下磁化強度的常微分方程：

17、

18、其中，參數α和k可通過磁滯回線中矯頑點(h＝he，m＝0)和剩磁點(h＝0，m＝mr)確定。在δ＝+1的矯頑點處，mja模型方程(4)為：

19、

20、同樣地，在δ＝-1的剩磁點，(4)可以寫成：

21、

22、通過組合式(5)和(6)可以確定參數α和k。求解常微分方程，mja和h之間的關系可以確定為：

23、

24、基于公式(7)，靈敏度為p(r)的接收線圈在r處感應電壓可以表示為：

25、

26、其中，m(r,t)＝c(r)·mja(r,t)，c(r)為位置r處的粒子濃度。定義s(r,t)為描述粒子濃度和電壓之間的映射的積分核：

27、

28、則公式(8)可以表示為：

29、u(t)＝∫ωc(r)·s(r,t)dr??(10)

30、為了進行數值計算，將公式(10)的時間與空間進行離散化：

31、

32、其中，s(rn,tk)的計算公式為：

33、

34、其中，m表示單個粒子的磁矩的模值。對公式(11)兩側在時域進行快速傅里葉變換(fft)得到：

35、

36、其中，將s(fk,rn)定義為系統矩陣，定義為mnps的實測電壓信號頻譜矢量。k＝1,2,3...,k為頻率分量個數，n＝1,2,3...,n為成像網格上的體素點個數。

37、進一步，步驟二具體方法步驟如下：

38、考慮磁納米粒子解空間的獨特性質及空間相鄰體素間磁化的內在關聯性，分別設計了兩種先驗信息。一是l1正則化先驗，旨在促進解的稀疏性，有效抑制噪聲干擾，提高重建圖像的質量；二是磁化權重的laplacian正則化先驗，該約束通過描述相鄰粒子空間中的磁化梯度關系，進一步增強了重建結果的局部連貫性和空間一致性。

39、首先，將公式(13)用矩陣的形式描述為：

40、

41、其中，為系統矩陣，為粒子濃度分布。為mnps的實測信號頻譜矢量。k為頻率分量個數，n為成像網格上的體素點個數。

42、第二，由于k遠小于n，公式(14)中對粒子濃度分布求解的問題可以轉換為求解不適定的逆問題，一般采用最小二乘法加正則化約束項求解，可以表示為：

43、

44、其中，為保證重構結果c盡可能接近真值的保真項。r(c)為正則化約束項。通過在反問題的保真項上引入附加約束，來克服反問題的病態性，從而獲得更穩定可靠的解。

45、在磁粒子成像技術中，感應電壓的產生主要依賴于零磁場區在視場(fov)內的動態移動，這一過程激發了磁納米粒子的磁化響應。值得注意的是，磁納米粒子在視場內的實際分布往往高度稀疏，即它們通常僅占據fov中極少數的體素點。這種高度的稀疏性特性對mpi圖像重建算法提出了特殊要求，需要設計能夠有效處理大規模數據稀疏性的方法，以確保從有限的磁化響應數據中精確重建出磁納米粒子的空間分布。因此，鑒于磁納米粒子在視場(fov)內分布的稀疏性，采用l1正則化項作為先驗約束條件。l1正則化項通過促進解的稀疏性，即減少非零元素的數量，有助于在重建過程中抑制噪聲的干擾，從而提高重建圖像的質量和準確性。對于公式(15)，r(c)可以表示如下：

46、

47、此外，必須指出的是，雖然促進稀疏性的正則化項能夠有效減少非零元素，從而保持解空間的稀疏性，但過度稀疏化可能會導致重建目標形態的損失，即關鍵細節和結構的模糊化。為了平衡稀疏性和重建質量，進一步提出結合空間先驗信息的策略，以確保重建圖像既具備形態學上的準確性，又不失細節豐富性。在mpi的空間編碼機制中，零磁場區域(ffr)內的磁納米粒子與鄰近粒子之間存在著顯著的時空相關性。特別是在中心體素點進行空間編碼的瞬間，其周邊體素點亦會產生微弱的磁化響應，這些響應在信號采集過程中被有效捕獲。基于這一高度相關性，采用一種基于拉普拉斯(laplacian)矩陣的磁化權重先驗約束項，該約束項旨在精確刻畫相鄰體素點之間的相互作用與依賴關系，從而在保證稀疏性的同時，最大限度地保留和恢復重建圖像的形態學特征和細節信息。

48、拉普拉斯矩陣的可以用l表示，其中l＝d-w，w是鄰接矩陣，表示每個體素點與其他每個點之間的強度關系。d為加權度矩陣，其中每個元素dii表示體素點i與其他點之間的總連接權重。通過這樣的定義，拉普拉斯矩陣能夠有效地刻畫圖像中體素點之間的局部連接性和差異性，為進一步利用空間先驗信息來約束和優化重建圖像的形態學和細節提供了數學基礎。首先，在mpi中一般認為如果施加的外加磁場使粒子的磁化強度達到飽和磁化強度的80％時，粒子達到磁化飽和，并將該外加磁場強度h的模值定義為飽和場強hsat。選擇場和驅動場疊加后，ffr處的粒子合磁場h(r,t)＝hd(t)+hs(r)＝0，則ffr周圍粒子rδ的合磁場可以表示為：

49、

50、其中，δt為ffr周圍坐標位置相對于ffr坐標位置處的位移矢量。當||h(rδ,t)||<hsat時，定義所有的構成的空間范圍為ffr的δ鄰域。對于鄰域內各體素點與ffr之間的關系權重，可以利用mja模型的導函數來衡量，因為由公式(12)知空間中每一個坐標位置r處產生的磁化響應與磁化強度m(r,t)對時間的變化率都成正比，其中因此，在ffr位置確定時，空間編碼的時刻tffr也是定值，mja(β||h(r,t)||)對h(r,t)的導數反映了在該時刻整個空間的磁化響應強度。鄰接矩陣w可以表示為：

51、

52、其中，i和j表示每一個體素點。δ(i)表示i點作為ffr點時該點的鄰域。δrij表示鄰域內點j到i的位移矢量。表示在h＝g·δrij時mja(β||h||)對h的導數值。表示在i點作為ffr點時mja(β||h||)對h的導數，它的值為β/3。對于δ(i)鄰域的確定，若在體素點i處外加磁場強度模值小于hsat，則該體素點就在ffr點i的鄰域δ(i)內，即通過求解0.8msat＝cmmja(β||h||)求解hsat。由于數值計算時體素點的坐標是離散的，因此要確定點j是否在i的鄰域內，只需要根據選擇場梯度g判斷||g·δrij||是否小于hsat即可。對于一個線性度良好的選擇場，各體素點磁化空間鄰域內的體素點的個數和分布是相同的，因此對于每個體素點鄰域內相對權重的計算可以簡化為只計算空間中一個體素點的鄰域內各點權重，再根據空間相對位置將該點鄰域的所有權重匹配到其他到體素點的磁化空間鄰域點上，從而完成了鄰接矩陣w的構建。對于加權度矩陣d可以表示為：

53、

54、根據公式(18)和公式(19)，拉普拉斯矩陣l的元素lij可以表示為：

55、lij＝dij-wij??(20)

56、根據d和w的定義，lij的值可以分三種情況：

57、1)當i＝j時：

58、

59、注意到wii＝0，因此，

60、

61、為了清晰起見，將求和中的索引從j改為了k，以避免與外部的j混淆。但在實際應用中，它們是相同的。

62、2)當i≠j且i和j相鄰時(j∈δ(i))：

63、lij＝dij-wij＝0-wij＝-wij??(23)

64、因此，

65、

66、3)當i≠j且i和j不相鄰時

67、lij＝dij-wij＝0-0＝0??(25)

68、綜合以上情況，拉普拉斯矩陣l的表達式為：

69、

70、對于n個單獨的體素點，l是一個大小為n×n關于主對角線對稱的對稱矩陣。至此，對于公式(15)，r(c)再加入拉普拉斯先驗約束項后可以進一步表示為：

71、

72、進一步，所述的步驟三包括以下步驟：

73、在步驟二所確立的先驗信息基礎上，構建了一個雙先驗約束的目標函數。此目標函數的核心作用需要平衡數據保真度(即確保重建結果與原始測量數據的一致性)與解的稀疏性、平滑性之間的微妙關系。因此，基于公式(27)的稀疏先驗約束項和空間先驗約束項，對于公式(15)可以寫為如下形式：

74、

75、其中，是濃度向量c的最優近似解，λ1和λ2是正則化參數，||c||1是l1正則化稀疏約束項，是拉普拉斯空間約束項。該目標函數根據mpi的成像特點引入稀疏約束以促進解的稀疏性和提高重建圖像的抗噪聲性能；引入基于磁化權重的laplacian正則項以保留像素點鄰域信息和恢復圖像的形態學信息。通過綜合這些正則化項來平衡重建圖像的稀疏性和形態學從而實現成像質量和空間分辨率的提高。

76、進一步，所述的步驟四包括以下步驟：

77、為了對目標函數公式(28)進行高效迭代求解，基于非凸優化和迭代逼近理論，本發明構建了深度展開模型，以達到對磁納米粒子濃度分布的高保真度重建。首先，進行梯度計算，對于公式(28)的第一項可以通過鏈式法則對c求導。展開該項為：

78、

79、對c求導的結果可以表示為：

80、

81、第二項α||c||1，其偏導數是非光滑的。l1范數的次梯度可以表示為：

82、

83、其中，sign(c)是符號函數，表示對每個元素的導數分別為+1(如果該元素為正)，-1(如果該元素為負)，以及0(如果該元素為0)。

84、第三項展開該項為：

85、ω(lc)t(lc)＝ωctltlc??(32)

86、對c求導的結果可以表示為：

87、

88、將各項偏導數結合起來，總的梯度公式可以表示為：

89、

90、基于非凸優化和迭代逼近理論，經過多次迭代可以找到迭代過程中更新的粒子信息的最優值。第k+1次迭代可以表示為：

91、c(k+1)＝c(k)-[st(sc(k)-u)+α·sign(c(k))+2ωltlc(k)]??(35)；

92、接下來，整個深度展開重建框架的核心構建是基于迭代公式(35)，它將迭代求解過程轉化為深度網絡的層級結構。每一次迭代都被映射為深度網絡中的一層，通過逐層迭代逼近指導著整個重建過程的推進。擴展后的深度神經網絡模型可以形式化地表示為：

93、c(k+1)＝ψθ(c(k)；u,s,l)??(36)

94、其中c(k)表示第k次迭代的濃度估計，而θ是所有迭代步驟中共享的模型參數集。此模型將拉普拉斯矩陣l和復數形式的測量電壓u及系統矩陣s作為輸入信息，粒子濃度c則在初始階段被全零初始化，旨在為后續的迭代過程提供一個基準起點。

95、模型的關鍵組成部分如下：

96、梯度下降優化策略：在每一迭代步驟中，采用梯度下降法來更新濃度估計c(k)。具體而言，通過計算保真項的負梯度方向，引導每次迭代的結果向更逼近真實解的方向移動。這一過程確保了迭代過程的收斂性和重建結果的準確性。

97、復數結果處理機制：鑒于公式(34)產生的結果為復數形式，采用通道和維度轉換層。這一層負責將每次迭代產生的一維復數輸出轉換為二維數據格式，并巧妙地分離實部和虛部，分別分配獨立的通道進行處理。這樣的設計不僅簡化了復數運算的復雜度，還提高了模型處理復數數據的能力。

98、特征提取與增強：為了在保證重建質量的同時有效抑制噪聲和冗余信息，引入了由編碼器、解碼器和和空間注意模塊構成的近端投影層。編解碼器作為特征提取器，能夠有效地捕捉和重構圖像中的關鍵特征表示。而空間注意模塊則進一步增強了模型對核心特征信息的敏感性，通過動態調整不同空間位置的權重，使模型能夠聚焦于對重建結果最為關鍵的區域，從而在抑制無關信息的同時提升重建精度。

99、模型的計算流程如下：

100、首先，將迭代優化過程顯式地表示為深度神經網絡的結構，從而利用梯度下降的基本思想對公式(35)進行優化。在第k次迭代中，計算保真度項相對于當前解項(記為c(k)的梯度其中表示保真度損失函數。隨后，使解項沿保真度項的負梯度方向移動一個小的步長η(學習率)，可以表示為：

101、

102、在這個過程中，迭代產生的結果用復數表示。其次，在通道和維度轉換層，將一維數據轉換為二維數據，并將實部和虛部分別提取到兩個不同的通道中。為后續的網絡層提供了更直觀、更易于處理的數據格式。轉換過程可以表示為：

103、

104、第二，為了在保證重建質量的同時有效抑制噪聲和冗余信息，近端投影層集成了編碼器、解碼器和空間注意模塊。編碼器部分采用批歸一化層來加速收斂和穩定梯度，并通過堆疊的卷積層和relu激活函數結構提高特征表示的復雜性和深度。卷積核的大小為3，輸出通道的數量依次為2,4,4,8,8,16，假設輸入的特征為編碼過程可以表示為：

105、

106、其中，hi是第i層卷積后的特征表示。接下來，引入空間注意模塊，為每個位置的粒子特征增加權重，進一步增強特征表征。

107、g6＝fsam(h6)??(40)

108、其中，fsam(·)表示空間注意力計算。然后，通過與編碼器結構相同的六層卷積操作解碼器逐步減少高級增強特征表示，設置反向映射通道的數量依次為16,8,8,4,4,2，可以表示為：

109、gi＝relu(conv(gi+1)),i＝5,4,…,0??(41)

110、其中，gi是第i層卷積后的特征表示。第四，為了進一步提高預測的準確性，平衡預測值向目標逼近的趨勢，通過反向通道和維度轉換層將解碼器的輸出轉換為復數形式，并根據迭代公式(35)使用殘差連接策略更新預測值。這一步有助于網絡在迭代過程中保留對初始值的記憶，并逐步提高預測的準確性。其中迭代次數預設為k＝5。最后，在整個迭代完成后，通過bn層去除內部協變量偏差來穩定和加速訓練過程，并通過兩次卷積運算，進一步提高網絡的非線性表示能力使預測結果近似于目標結果。

111、此外補充說明，涉及到的空間注意模塊fsam(·)接收一個特征圖x作為輸入，其尺寸為c×h×w，其中c是通道數，h和w分別是特征圖的高度和寬度，該特征圖是的輸出，輸出的結果作為的輸入。模塊通過計算空間上的權重圖，來重新加權輸入特征圖的不同位置，從而增強對重要區域的關注。首先，使用全局平均池化和全局最大池化將每個通道的空間信息壓縮為一個全局標量。

112、

113、其中，favg和fmax本別表示全局平均池化和全局最大池化。接著，通過卷積層對這兩個向量進行變換，以提取更復雜的特征。變換后的向量記為favg和favg。

114、

115、其中，wavg和wmax是權重矩陣，bavg和bmax是偏置項，σ是relu非線性激活函數。然后，將變換后的向量進行元素級的拼接，并通過一個額外的卷積和一個sigmoid激活函數，生成空間注意力權重圖m，其尺寸為1×h×w。

116、

117、最后，將注意力權重圖m與原始輸入特征圖x進行逐元素相乘，得到加權后的特征圖x′。

118、x′(c,i,j)＝x(c,i,j)×m(i,j)??(45)；

119、進一步，所述的步驟五包括以下步驟：

120、采用均方誤差(mse,mean?squarederror)和結構相似性指數(ssim,structuralsimilarityindex)的加權和作為損失函數用于網絡模型訓練，總的損失函數可以表示為：

121、loss＝λ·mse+(1-λ)·(1-ssim)??(46)

122、其中，λ＝0.9是一個權重參數，用于平衡mse和ssim之間的相對重要性。mse用于評估重建圖像的質量。ssim用于評估重建圖像與原始圖像之間的相似性。mse可以表示為：

123、

124、其中，y和分別為原始和網絡預測的mpi圖像，m為圖像個數。ssim可以表示為：

125、

126、其中，y和分別為原始和網絡預測的mpi圖像。μy和分別是圖像y和的均值。和分別是圖像y和的方差。是圖像y和的協方差。c1和c2是用于穩定性的常數，對于歸一化圖像設置為0.01和0.03。然后，利用訓練的深度神經網絡模型，對測試數據集進行測試，以驗證模型在實際應用中的性能。同時，展示了預測結果和與真實圖像之間的誤差圖。

127、與現有技術相比，本發明具有以下優點：

128、第一，本發明引入改進的jiles-atherton磁滯模型來代替朗之萬函數，以更準確地模擬spios在mpi系統中的磁化行為，從而提高測量電壓和系統矩陣的模擬準確度。

129、第二，本發明將傳統的迭代算法與深度學習相結合，通過訓練隱函數的方式優化算法性能。深度學習參數訓練減輕手動調參的難度，以自適應不同的圖像重建任務和需求。

130、第三，本發明考慮了mpi磁化空間內鄰近體素點之間的空間相關性，這是傳統迭代算法所忽視的。通過引入基于磁化權重的laplacian正則化項，能夠更準確地描述和模擬磁化空間的相關性，從而在一定程度上提高圖像重建的準確性和魯棒性。