Featured image of post 從陰影與投影能重構看不見的內部嗎?:反問題、拉東變換與CT斷層成像的數學

從陰影與投影能重構看不見的內部嗎?:反問題、拉東變換與CT斷層成像的數學

如何僅憑從外部穿透物體的X射線投影數據,重構出人體或物體內部精確的三維解剖結構?本文深入剖析反問題(Inverse Problems)的本質、阿達馬適定性條件、拉東變換(Radon Transform)的積分幾何理論、傅立葉切片定理與濾波反投影(FBP),以及應對不適定性災難的吉洪諾夫正規化與壓縮感知前沿。

導言:透過陰影窺探不可見的深處

在日常經驗中,因果邏輯的推進往往是單向且自然的:如果知道了物體的形狀和光源的位置,投射在牆面上的陰影輪廓便能通過幾何光學極其容易地計算出來;如果知道了地層的物質密度分佈與彈性參數,地震波從震源傳播到地面各個地震儀的走時與波形便能根據波動方程精確模擬。這種從原因(內在系統參數與幾何構造)推導結果(外部觀測數據與物理場響應)的過程,被稱為正問題(Forward Problem)。正問題通常對應著大自然本身的物理運行規律,在數學上往往具有確定、穩定且唯一的解。

然而,在現代科學、醫學診斷和工程勘探的核心場景中,我們面臨的往往是截然相反的逆向挑戰:

  • 醫學電腦斷層掃描(CT / Computed Tomography):醫生無法切開患者的顱骨或胸腔,只能從外部不同角度發射X射線並測量衰減後的透射強度。我們能否根據這些一維或二維的投影「陰影」,精確重構出大腦軟組織、微小腫瘤或血管壁的內部三維密度分佈?
  • 地球物理勘探:地質學家在地面引爆炸藥或震源車,記錄地面各處的地震波反射信號。我們能否據此反推地下數千米處石油構造、斷層滑移面或地幔對流邊界?
  • 天文干涉成像:事件視界望遠鏡(EHT)散佈在地球各處的射電天線僅採集到了極其稀疏的空間相干度數據。我們能否逆向合成出超大質量黑洞事件視界周圍吸積盤的真實圖像?
  • 無損檢測與材料科學:工業超聲與中子斷層成像如何在不破壞航空發動機單晶渦輪葉片的前提下,發現其微米級內部微裂紋?

這些從觀測結果(數據)逆向推導物理原因(系統內部狀態)的數學任務,統稱為反問題(Inverse Problems)。反問題不僅是應用數學中最具智力挑戰和美感的領域之一,更是現代無損成像技術的絕對理論支柱。

本文將帶領讀者深入反問題與現代成像的數學殿堂:從阿達馬(Jacques Hadamard)定義的「不適定性災難」,到約翰·拉東(Johann Radon)於1917年奠基的積分幾何與拉東變換;從傅立葉切片定理(Fourier Slice Theorem)構建的對偶橋樑,到濾波反投影演算法(FBP)的優雅推導;最終探索吉洪諾夫正規化(Tikhonov Regularization)、全變差正規化(Total Variation)以及壓縮感知(Compressed Sensing)如何馴服高維噪聲與病態矩陣,將不可能變為現實。

  flowchart LR
    subgraph Forward ["正問題 (Forward Problem)"]
        direction TB
        A["物理本質 f(x, y)<br/>(內部組織密度)"] -->|物理定律 / 算子 A| B["觀測數據 g(s, θ)<br/>(X射線投影陰影)"]
    end
    subgraph Inverse ["反問題 (Inverse Problem)"]
        direction TB
        C["外部觀測數據 g(s, θ)<br/>(含噪聲觀測)"] -->|逆算子 A⁻¹ / 正規化| D["重構內部結構 f*(x, y)<br/>(高精度斷面圖像)"]
    end

1. 什麼是反問題?:阿達馬適定性與不適定性的災難

要理解反問題的深刻性,首先必須明晰為什麼求解反問題比求解正問題困難千百倍。

1.1 算子方程的形式化表述

$$g = \mathcal{A}(f)$$

其中:

  • $f \in \mathcal{X}$ 表示未知的系統內部物理參數分佈(例如人體內部的X射線線性衰減係數分佈 $f(x, y) \ge 0$);
  • $g \in \mathcal{Y}$ 表示觀測儀器所收集到的外部測量數據(例如探測器陣列測得的光子衰減積分投影數據);
  • $\mathcal{A}: \mathcal{X} \to \mathcal{Y}$ 為正向物理算子。對於CT成像、熱傳導或線性波動而言,$\mathcal{A}$ 通常是一個連續的線性積分算子。
$$\mathcal{A} f = g$$

1.2 阿達馬適定性準則(Hadamard Well-Posedness)

法國數學家雅克·阿達馬在1902年系統研究偏微分方程邊值問題時提出:一個數學物理問題被稱為適定問題(Well-Posed Problem),當且僅當滿足以下三條基本準則:

  1. 解的存在性(Existence):對於任意可能的數據 $g \in \mathcal{Y}$,方程 $\mathcal{A} f = g$ 至少存在一個解 $f \in \mathcal{X}$。
  2. 解的唯一性(Uniqueness):方程 $\mathcal{A} f = g$ 的解至多只有一個。若 $\mathcal{A} f_1 = \mathcal{A} f_2$,則必然有 $f_1 = f_2$(即算子 $\mathcal{A}$ 必須是單射,其核空間 $\ker(\mathcal{A}) = \{0\}$)。
  3. 解的連續依賴性 / 穩定性(Stability):解 $f$ 連續依賴於輸入數據 $g$。具體而言,逆映射 $\mathcal{A}^{-1}: \mathcal{Y} \to \mathcal{X}$ 是連續算子。如果測量數據產生微小的擾動 $\delta g$(如儀器熱噪聲、光子散粒噪聲),重構解的誤差 $\delta f = \mathcal{A}^{-1}(g + \delta g) - \mathcal{A}^{-1}(g)$ 應當能夠被常數控制,即: $$\|\delta f\|_{\mathcal{X}} \le C \|\delta g\|_{\mathcal{Y}}$$

如果上述三個條件中任何一個不滿足,該數學問題就被稱為不適定問題(Ill-Posed Problem)。

1.3 反問題的不適定性根源:緊算子與高頻災難

令人驚異的是,現實物理世界中的幾乎所有反問題,天生都是嚴重不適定(Severely Ill-Posed)的,其根本原因在於物理正向過程的平滑效應(Smoothing Effect)。

$$g(s) = (\mathcal{A} f)(s) = \int_a^b K(s, t) f(t) \, dt$$

其中核函數 $K(s, t)$ 通常是光滑的高斯核、泊松核或波動格林函數。

根據積分算子的性質,積分操作本質上是一種低通加權平均。無論未知的原信號 $f(t)$ 包含多麼劇烈的高頻振盪、尖銳邊緣或階躍斷層,經過連續核 $K(s, t)$ 的平滑積分之後,生成的數據 $g(s)$ 都會變得異常光滑。

$$\sigma_1 \ge \sigma_2 \ge \dots \ge \sigma_n \to 0 \quad (n \to \infty)$$$$\mathcal{A} f = \sum_{n=1}^\infty \sigma_n \langle f, v_n \rangle u_n \implies f = \mathcal{A}^{-1} g = \sum_{n=1}^\infty \frac{1}{\sigma_n} \langle g, u_n \rangle v_n$$$$\delta f = \frac{\epsilon}{\sigma_n} v_n$$

當模態階數 $n$ 增大時,$\sigma_n$ 指數級或多項式級衰減至 0,導致放大的噪聲項 $\frac{\epsilon}{\sigma_n} \to \infty$ 徹底吞噬整個解!沒有任何直接求逆算子能夠在不加干預的情況下從含噪數據中恢復真實物理量。


2. 積分幾何的先驅:約翰·拉東與拉東變換

在X射線衍射和醫學CT問世前半個世紀,奧地利數學家約翰·拉東(Johann Radon)於1917年在純粹微分幾何與數學分析的研究中,提出了一個極為優美的幾何問題:

「如果我們知道了定義在二維平面上的未知光滑函數沿著所有可能直線路徑的線積分,我們能否唯一且精確地確定原函數本身?」

拉東不僅證明了其唯一性,更推導出了嚴格的解析逆變換公式。這一純數學創舉在半個世紀後成為了價值百億美元的現代醫學CT工業的數學地基。

2.1 平面直線的法線參數化

$$L(s, \theta) = \left\{ (x, y) \in \mathbb{R}^2 \;\middle|\; x \cos\theta + y \sin\theta = s \right\}$$$$\mathbf{x}(t) = s \boldsymbol{\omega} + t \boldsymbol{\omega}^\perp = (s \cos\theta - t \sin\theta, \; s \sin\theta + t \cos\theta)$$

其中 $t \in (-\infty, \infty)$ 為直線上的局部線坐標。

2.2 拉東變換(Radon Transform)的數學定義

$$(\mathcal{R}f)(s, \theta) = \int_{-\infty}^\infty f(s\cos\theta - t\sin\theta, \; s\sin\theta + t\cos\theta) \, dt$$$$(\mathcal{R}f)(s, \theta) = \iint_{\mathbb{R}^2} f(x, y) \, \delta(x \cos\theta + y \sin\theta - s) \, dx dy$$$$s(\theta) = x_0 \cos\theta + y_0 \sin\theta = r_0 \cos(\theta - \phi_0)$$

3. 醫學CT的物理本質:比爾-朗伯定律與X射線投影

為什麼拉東變換能直接對應現代醫院裡的CT機工作原理?這一切源於經典的電磁輻射吸收物理學。

3.1 比爾-朗伯定律(Beer-Lambert Law)

$$\frac{dI(t)}{dt} = - f(\mathbf{x}(t)) \, I(t)$$$$\int_{I_0}^{I_{\text{det}}} \frac{dI}{I} = - \int_{L} f(\mathbf{x}) \, dt$$$$I_{\text{det}} = I_0 \exp\left( - \int_L f(\mathbf{x}) \, dt \right)$$

其中 $I_0$ 為X射線管發出的初始光強。

3.2 對數變換:從非線性測量到拉東投影

$$g(s, \theta) = \ln\left( \frac{I_0}{I_{\text{det}}} \right) = \int_L f(\mathbf{x}) \, dt = (\mathcal{R}f)(s, \theta)$$

因此,醫學CT圖像重構問題,在數學上等價於:已知全角度的拉東變換正弦圖數據 $g(s, \theta) = (\mathcal{R}f)(s, \theta)$,求解未知的人體內部衰減係數分佈 $f(x, y)$。


4. 傅立葉切片定理(Fourier Slice Theorem):通往解析求逆的鑰匙

如何從一維投影積分 $g(s, \theta)$ 恢復出二維圖像 $f(x, y)$?連接投影空間與物理圖像空間的極其關鍵的數學橋樑,就是傅立葉切片定理(Fourier Slice Theorem / 投影切片定理)。

4.1 定理的數學陳述

定理:未知二維函數 $f(x, y)$ 的二維連續傅立葉變換,在頻域極坐標系下沿角度 $\theta$ 的一維徑向切片,嚴格等於該函數在同角度 $\theta$ 下的一維投影拉東變換 $g(s, \theta)$ 關於探測器位移坐標 $s$ 的一維傅立葉變換!

$$\left( \mathcal{F}_1 [\mathcal{R}f(\cdot, \theta)] \right)(\omega) = \left( \mathcal{F}_2 f \right)(\omega \cos\theta, \; \omega \sin\theta)$$

4.2 嚴密數學推導

$$G(\omega, \theta) = \int_{-\infty}^\infty g(s, \theta) \, e^{-i 2\pi \omega s} \, ds$$$$G(\omega, \theta) = \int_{-\infty}^\infty \left[ \iint_{\mathbb{R}^2} f(x, y) \, \delta(x \cos\theta + y \sin\theta - s) \, dx dy \right] e^{-i 2\pi \omega s} \, ds$$$$G(\omega, \theta) = \iint_{\mathbb{R}^2} f(x, y) \left[ \int_{-\infty}^\infty \delta(x \cos\theta + y \sin\theta - s) \, e^{-i 2\pi \omega s} \, ds \right] dx dy$$$$G(\omega, \theta) = \iint_{\mathbb{R}^2} f(x, y) \, e^{-i 2\pi \omega (x \cos\theta + y \sin\theta)} \, dx dy$$$$G(\omega, \theta) = \iint_{\mathbb{R}^2} f(x, y) \, e^{-i 2\pi (u x + v y)} \, dx dy = \hat{f}(u, v) = \hat{f}(\omega \cos\theta, \omega \sin\theta)$$

證畢!


5. 濾波反投影演算法(FBP):工業CT的經典皇冠

擁有傅立葉切片定理後,推導出更加穩健、高效且可即時串流計算的演算法——濾波反投影演算法(Filtered Backprojection / FBP)。

5.1 二維逆傅立葉變換的極坐標分解

$$f(x, y) = \iint_{\mathbb{R}^2} \hat{f}(u, v) \, e^{i 2\pi (u x + v y)} \, du dv$$$$u = \omega \cos\theta, \quad v = \omega \sin\theta$$$$du dv = |\omega| \, d\omega d\theta$$$$f(x, y) = \int_0^\pi d\theta \int_{-\infty}^\infty \hat{f}(\omega \cos\theta, \omega \sin\theta) \, e^{i 2\pi \omega (x \cos\theta + y \sin\theta)} \, |\omega| \, d\omega$$$$f(x, y) = \int_0^\pi \left[ \int_{-\infty}^\infty G(\omega, \theta) \cdot |\omega| \cdot e^{i 2\pi \omega s} \, d\omega \right]_{\!s = x\cos\theta + y\sin\theta} d\theta$$

5.2 濾波與反投影的兩步解耦

整個重構過程被優美地解耦為兩步操作:

第一步:斜坡濾波(Ramp Filtering)

$$\tilde{g}(s, \theta) = \mathcal{F}_1^{-1} \left\{ G(\omega, \theta) \cdot |\omega| \right\} = g(s, \theta) * h(s)$$

其中濾波器頻域傳遞函數為 $H(\omega) = |\omega|$,被稱為Ram-Lak(斜坡)濾波器。它精準地抵消了 $1/r$ 的幾何反投影模糊!

第二步:反投影累加(Backprojection)

$$f(x, y) = \int_0^\pi \tilde{g}(x \cos\theta + y \sin\theta, \; \theta) \, d\theta = \mathcal{B} \tilde{g}$$

6. 不適定性的馴服:吉洪諾夫正規化與變分方法

當投影角度不足、噪聲劇烈或劑量極低時,解析算法失效,必須轉向現代正規化技術。

6.1 離散代數重構模型

$$\mathbf{A} \mathbf{f} + \boldsymbol{\eta} = \mathbf{g}$$

其中矩陣元 $A_{ij}$ 表示第 $i$ 條射線穿過第 $j$ 個像素的相交長度。

6.2 吉洪諾夫正規化(Tikhonov Regularization)

$$\mathbf{f}_{\alpha} = \arg\min_{\mathbf{f}} \left\{ \|\mathbf{A} \mathbf{f} - \mathbf{g}\|_2^2 + \alpha \|\mathbf{L} \mathbf{f}\|_2^2 \right\}$$

其中第一項為數據保真項,第二項為正規化懲罰項,$\alpha > 0$ 為調節權衡的正規化參數。

$$\mathbf{f}_\alpha = (\mathbf{A}^T \mathbf{A} + \alpha \mathbf{I})^{-1} \mathbf{A}^T \mathbf{g}$$$$\mathbf{f}_\alpha = \sum_{i=1}^r \left( \frac{\sigma_i^2}{\sigma_i^2 + \alpha} \right) \cdot \frac{\mathbf{u}_i^T \mathbf{g}}{\sigma_i} \mathbf{v}_i$$

濾波因子 $\frac{\sigma_i^2}{\sigma_i^2 + \alpha}$ 在高頻噪聲模態($\sigma_i \to 0$)時強行將增益截斷,避免了數值崩潰!


7. 現代前沿:全變差(TV)、壓縮感知與深度學習反演

7.1 全變差正規化(Total Variation)與邊緣保持

$$\min_{\mathbf{f}} \left\{ \frac{1}{2} \|\mathbf{A}\mathbf{f} - \mathbf{g}\|_2^2 + \lambda \|\nabla \mathbf{f}\|_1 \right\}$$

在消除噪聲的同時,完美保留了腫瘤與器官邊緣的銳利跳變。

7.2 壓縮感知理論(Compressed Sensing)

陶哲軒與坎德斯證明:只要原圖像在某種基底下稀疏,即使僅採樣極少數稀疏角度,也能無失真重構,使人體承受的X射線輻射劑量劇降80%以上!

7.3 深度學習與物理知情模型(Deep Unfolding)

將經典近端梯度或ADMM演算法展開為深度網絡,把物理正向算子 $\mathcal{A}$ 作為嚴格先驗層內嵌,兼具深度學習的高速度與物理反演的嚴格可靠性。


結語:在不可見的投影中構建理性之光

從柏拉圖洞穴寓言中的陰影,到現代醫學以微米精度切片人體生命構造;反問題的數學歷程展现了人類理性的偉大光輝:即使肉眼無法直接透視內部,憑藉嚴謹的積分幾何、泛函分析與數值演算法,我們依然能夠從外部有限的投影陰影中,毫釐不差地重構出深藏內部的真實全貌。

comments powered by Disqus