以下我們討論 為何 線性非時變 (Linear Time-Invariant, LTI) 系統 輸入與輸出關係 由 所謂的 convolution 表示。為了避免過多繁雜的數學，以下僅討論離散時間的情況。首先我們需要一些定義的幫助：
========================
Definition: Unit Impulse Function in Discrete-Time (Kronecker Delta)
我們說函數 $\delta: \mathbb{N} \to \{0,1\}$ 為 unit impulse function in discrete-time time 若 $\delta$ 滿足
\[\delta \left[ n \right] = \left\{ \begin{gathered}
1,\;\;\;\;n = 0 \hfill \\
0,\;\;\;\; o.w. \hfill \\
\end{gathered} \right.
\]========================
========================
Definition: Impulse Response
給定任意系統配備輸入 $x[n]$ 與輸出 $y[n]$ 關係為 $y[n] = T\{x[n]\}$其中 $T$ 視為 operator (定義在某函數空間)，若輸入為 $x[n]=\delta[n]$ 則 輸出
\[
h[n] := y[n] = T\{\delta[n]\}
\]稱為系統 $T$ 的 脈衝響應 (impulse response)
========================
========================
FACT: 任意離散訊號 $x[n]$ 可由 $\delta$ 做組合疊加，亦即
\[
x[n] = \sum_{k=-\infty}^{\infty}x[i] \delta[n-k]
\]========================
Proof: 證明顯然，在此不做贅述。$\square$
========================
Definition: Linear System
給定系統 $T$ 滿足以下輸入與輸出關係： $y_1[n]=T\{x_1[n]\}$ 且 $y_2[n] = T\{x_2[n]\}$。現在定義 $x[n] =ax_1[n] + bx_2[n]$ 其中 $a,b \in \mathbb{R}$則我們說系統 $T$ 為 linear 若下列條件成立
\[
y[n] = T\{x[n]\}
\]且 $y[n] = a y_1[n] + b y_2[n]$
========================
Remarks:
上述 系統 $T$ 其實可看按作 泛函分析中的 operator。故線性系統 就是 泛函分析中的 線性算子 (linear operator)。
========================
Definition: Time-Invariant System
給定系統 $T$ 滿足以下輸入與輸出關係： $y[n]=T\{x[n]\}$。現在定義 $x[n] :=x[n-d]$ 其中 $d \in \mathbb{N}$ 則我們說 $T$ 為 time-invariant 若下列條件成立
\[
y[n-d] = T\{x[n-d]\}
\]========================
有了以上定義的幫助，我們現在有辦法給出以下為何 LTI 系統的輸入輸出關係確實為 convolution。
========================
Theorem: LTI input-output Relationship is Governed by Convolution
給定 線性非時變 (LTI)系統 $L$ 滿足
\[
y[n] = L\{x[n]\}
\]其中 $L$ 為 linear operator，則
\[
y[n] = \sum_k h[k] x[n-k]
\]========================
Proof:
由於 任意 輸入 $x[n]$ 皆可被表為 脈衝函數的 組合：
\[
x[n] = \sum_{k = -\infty}^\infty x[k] \delta[n-k]
\]現在將此訊號作為LTI系統 $L$ 的輸入，則由於 $L$ 為 linear operator 我們有
\begin{align*}
y[n] &= L\left\{ {x\left[ n \right]} \right\} \\
&= L\left\{ {\sum\limits_{k = - \infty }^\infty x [k]\delta [n - k]} \right\} \hfill \\
&= \sum\limits_{k = - \infty }^\infty x[k] L\left\{ \delta [n - k] \right\} \;\;\;\; (*)
\end{align*} 接著因為 $L$ 為 time-invariant 我們可進一步改寫
\[L\left\{ {\delta [n - k]} \right\} = h\left[ {n - k} \right] \;\;\;\; (**)
\]故由 $(*)$ 與 $(**)$ 可得
\[y[n] = \sum\limits_{k = - \infty }^\infty {x[k]h\left[ {n - k} \right]} \]即為所求。$\square$
Comments:
上述 \[
y[n] = \sum\limits_{k = - \infty }^\infty {x[k]h\left[ {n - k} \right]}
\]一般又稱作 $x[n]$ 與 $h[n]$ convolution sum。記作
\[
y[n] = x[n] * h[n]
\]
7/07/2017
[控制理論] 非線性模型預測控制 (0) - 引論
此文我們針對 非線性模型預測控制 (Nonlinear Model Predictive Control, NMPC) 做簡單的介紹： 是一套針對受控廠為 非線性系統 所發展的 回授控制理論，其控制律透過不斷求解最佳化問題而得。以下我們將其簡稱 NMPC。
基本想法：
假設你是一個圍棋高手正在與旗鼓相當的對手對弈，那麼你可能在心中盤算好幾步可能的走法 並同時 試圖來 "預測" 對手的路數，但輪到自己下子的時候，我們只能選取剛才在心中盤算出的所有走法中 " 最佳" 的那一步，並且只動那一步旗，動了之後都必須重新盤算上述過程。這個精神大概與 模型預測控制 本質幾乎相同。就是逐步最佳化，慢慢朝向目標前進。以下我們會逐步將此概念規範化，讀者可以讀完之後再回頭瞧瞧這個基本想法，也許會發現異曲同工之妙。
對 $n=0,1,2,...$，令 $x(n)$ 為當前系統狀態 且 $x^{ref}(n)$ 為 給定的參考目標軌跡。
模型預測控制 的 主要目的：
與一般控制問題相仿，模型預測控制的主要目的不外乎以下兩大類問題：
1. 鎮定(stabilizing)問題。使預測的狀態軌跡 盡可能趨近 零點。
2. 追蹤(tracking)問題。使 預測的狀態軌跡 盡可能跟隨給定的 參考目標軌跡 。
Comments:
1. 若稍有控制背景的讀者不難發現鎮定問題其實是追蹤問題的特例，更進一步地說，所謂追蹤控制問題是指：決定一組控制輸入 $u(n)$ 使得 $x(n)$ 能盡可能緊跟給定的參考狀態軌跡 $x^{ref}(n)$。換而言之，若當前狀態 $x(n)$ 與 $x^{ref}(n)$ 所差甚遠，則我們將盡可能控制此系統軌跡趨近 $x^{ref}(n)$：若當前狀態 $x(n) = x^{ref}(n)$ 則我們將盡可能控制此系統"維 持"在該狀態。
2. 關於鎮定問題提及的對零點追蹤 或者 追蹤問題的參考目標軌跡，都必須滿足隱藏的前提： 零點 或者 參考目標軌跡 必須是系統的 平衡點 (equilibrium point)。否則該系統無法進行鎮定 or 追蹤。
以下我們將介紹一類簡化的 NMPC 問題，並藉此展示此控制方法的一些基本精神。
一類簡化的非線性模型預測控制問題：
令 $x(n) \in X := \mathbb{R}^d$ 且 $u(n) \in U := \mathbb{R}^m$ 考慮 參考軌跡狀態(reference state trajectory)，記作 $x^{ref}(\cdot)$，寫作
\[
x^{ref}(n) :=x_* :=0,\;\;\;\;\; n\geq 0
\]上述 對零點狀態追蹤問題 退化為 鎮定問題。
以控制理論的基本想法，我們很自然想要對 當前狀態 $x(n)$ 與 目標 $x^{ref} = 0$ 之誤差 做回授控制，為了達成此目標，我們期望控制力 $u$ 必須具備以下形式：
$$
u(n) := \mu (x(n))
$$其中 $\mu $ 為 將狀態 $x \in X$ 映到控制輸出值之集合 $U$ 之映射。
系統模型動態(System Model Dynamics)：
模型預測控制的想法在於利用 系統的模型來進行未來狀態軌跡之預測，並且給出適當的目標函數後對其做最佳化求解控制力。更精確的說，系統動態 一般我們寫成下列方程
\[
x^+ = f(x,u) \;\;\;\;\; (*)
\]其中 $f: X \times U \to X$ 為已知且非線性 滿足 將狀態 $x$ 與控制輸出值 $u$ 映到下一個時刻的狀態 $x^+$ 。
預測模型動態(Predictive Model Dynamics)：
現在給定當前狀態 $x(n)$，對任意給定控制序列
\[
u(0),u(1),...,u(N-1)
\]滿足控制區間 $N \geq 2$，我們可以用此 控制序列 以及 給定的當前狀態 作為初始狀態來迭代 $(*)$，如此可得 預測狀態軌跡 $x_u$ 滿足下式：
\[
x_u(0) = x(n); \;\;\;\; x_u(k+1) = f(x_u(k), u(k)),\;\;\;\; k=0,...,N-1 \;\;\;\; (**)
\]由於 $x(n)$ 被作為初始狀態來進行預測，故我們 $x_u(k)$ 實際上為狀態 $x(n+k)$ (對時刻 $t_{n+k}$ ) 的 預測狀態。因此我們得到在時刻 $t_n, t_{n+1},...,t_{n+N}$ 對系統狀態的預測 (與狀態輸入序列 $u(0),...,u(N-1)$有關)。
每時刻的最佳控制力求解：
現在我們使用最佳控制的想法來決定 $u(0),...,u(N-1)$ 使得預測狀態 $x_u$ 能盡可能靠近 $x^{ref} = 0$。為了達成此目標，我們可選定 一距離函數
\[
l(x_u(k),u(k)) : = ||x_u(k)||^2 + \lambda ||u(k)||^2
\]其中 $||\cdot||$ 為 Euclidean norm 且 $\lambda \geq 0$為控制力的權重。則現在我們可以給出最佳控制問題如下
\[\min J\left( {x\left( n \right),u\left( . \right)} \right): = \sum\limits_{k = 0}^{N - 1} {l\left( {{x_u}\left( k \right),u\left( k \right)} \right)} \]對 $u(0),...,u(N-1)$ 為admissible 且 $x_u$ 滿足 (**)。
讓我們假設上述最佳控制問題具有解，且此解由一組最佳控制數列 記作
\[
u^*(0),...,u^*(N-1)
\]亦即我們有
\[\mathop {\min }\limits_{u\left( 0 \right),...,u\left( {N - 1} \right)} J\left( {x\left( n \right),u\left( . \right)} \right) = \sum\limits_{k = 0}^{N - 1} {l\left( {{x_{{u^*}}}\left( k \right),{u^*}\left( k \right)} \right)} \]現在為了讓控制力 $u$ 具有我們需要的回授形式 $\mu(x(n))$ 我們令
\[
\mu(x(n) ) := u^*(0)
\]亦即我們只使用最佳控制序列的第一元素。接著我們重複上述過程：在時刻 $t_{n+1}$ 可量得 $x(n+1)$ 並依此再度執行最佳化求得 $\mu(x(n+1))$。同理，對於時刻 $t_{n+2}$ 可量得 $x_{n+2}$ 並依此執行最佳化得到 $\mu(x(n+2))$...
Comments:
1. 因為受控系統為非線性，一般而言為我們所得到的 回授控制律 $\mu(\cdot)$ 必須透過迭代最佳化演算 而得，沒有解析解，故計算複雜度會是模型預測控制的一大挑戰。
2. 就預測模型觀點而言，我們可看出 預測狀態軌跡 $x_u(k), \;\; k=0,1,2,...,N$ 提供了對原本系統在已知時刻 $t_n$ 對於 $t_n,...,t_{n+N}$ 的預測，以及在時刻 $t_{n+1}$ 對於 $t_{n+1},...,t_{n+N+1}$ 的預測，以及 在時刻 $t_{n+2}$ 對於 $t_{n+2},...,t_{n+N+2}$ 的預測。因此我們可看出預測區間是 "移動" 的，故模型預測控制一般又稱為 移動區間控制 (Moving Horizon Control, MHC) 或稱 (Receding Horizon Control, RHC)
3. 若上述 $f$ 為線性，則上述討論稱之為 線性模型預測控制(Linear Model Predictive Control, LMPC)或者簡稱 模型預測控制 (Model Predictive Control, MPC)
基本想法：
假設你是一個圍棋高手正在與旗鼓相當的對手對弈，那麼你可能在心中盤算好幾步可能的走法 並同時 試圖來 "預測" 對手的路數，但輪到自己下子的時候，我們只能選取剛才在心中盤算出的所有走法中 " 最佳" 的那一步，並且只動那一步旗，動了之後都必須重新盤算上述過程。這個精神大概與 模型預測控制 本質幾乎相同。就是逐步最佳化，慢慢朝向目標前進。以下我們會逐步將此概念規範化，讀者可以讀完之後再回頭瞧瞧這個基本想法，也許會發現異曲同工之妙。
對 $n=0,1,2,...$，令 $x(n)$ 為當前系統狀態 且 $x^{ref}(n)$ 為 給定的參考目標軌跡。
模型預測控制 的 主要目的：
與一般控制問題相仿，模型預測控制的主要目的不外乎以下兩大類問題：
1. 鎮定(stabilizing)問題。使預測的狀態軌跡 盡可能趨近 零點。
2. 追蹤(tracking)問題。使 預測的狀態軌跡 盡可能跟隨給定的 參考目標軌跡 。
Comments:
1. 若稍有控制背景的讀者不難發現鎮定問題其實是追蹤問題的特例，更進一步地說，所謂追蹤控制問題是指：決定一組控制輸入 $u(n)$ 使得 $x(n)$ 能盡可能緊跟給定的參考狀態軌跡 $x^{ref}(n)$。換而言之，若當前狀態 $x(n)$ 與 $x^{ref}(n)$ 所差甚遠，則我們將盡可能控制此系統軌跡趨近 $x^{ref}(n)$：若當前狀態 $x(n) = x^{ref}(n)$ 則我們將盡可能控制此系統"維 持"在該狀態。
2. 關於鎮定問題提及的對零點追蹤 或者 追蹤問題的參考目標軌跡，都必須滿足隱藏的前提： 零點 或者 參考目標軌跡 必須是系統的 平衡點 (equilibrium point)。否則該系統無法進行鎮定 or 追蹤。
以下我們將介紹一類簡化的 NMPC 問題，並藉此展示此控制方法的一些基本精神。
一類簡化的非線性模型預測控制問題：
令 $x(n) \in X := \mathbb{R}^d$ 且 $u(n) \in U := \mathbb{R}^m$ 考慮 參考軌跡狀態(reference state trajectory)，記作 $x^{ref}(\cdot)$，寫作
\[
x^{ref}(n) :=x_* :=0,\;\;\;\;\; n\geq 0
\]上述 對零點狀態追蹤問題 退化為 鎮定問題。
以控制理論的基本想法，我們很自然想要對 當前狀態 $x(n)$ 與 目標 $x^{ref} = 0$ 之誤差 做回授控制，為了達成此目標，我們期望控制力 $u$ 必須具備以下形式：
$$
u(n) := \mu (x(n))
$$其中 $\mu $ 為 將狀態 $x \in X$ 映到控制輸出值之集合 $U$ 之映射。
系統模型動態(System Model Dynamics)：
模型預測控制的想法在於利用 系統的模型來進行未來狀態軌跡之預測，並且給出適當的目標函數後對其做最佳化求解控制力。更精確的說，系統動態 一般我們寫成下列方程
\[
x^+ = f(x,u) \;\;\;\;\; (*)
\]其中 $f: X \times U \to X$ 為已知且非線性 滿足 將狀態 $x$ 與控制輸出值 $u$ 映到下一個時刻的狀態 $x^+$ 。
預測模型動態(Predictive Model Dynamics)：
現在給定當前狀態 $x(n)$，對任意給定控制序列
\[
u(0),u(1),...,u(N-1)
\]滿足控制區間 $N \geq 2$，我們可以用此 控制序列 以及 給定的當前狀態 作為初始狀態來迭代 $(*)$，如此可得 預測狀態軌跡 $x_u$ 滿足下式：
\[
x_u(0) = x(n); \;\;\;\; x_u(k+1) = f(x_u(k), u(k)),\;\;\;\; k=0,...,N-1 \;\;\;\; (**)
\]由於 $x(n)$ 被作為初始狀態來進行預測，故我們 $x_u(k)$ 實際上為狀態 $x(n+k)$ (對時刻 $t_{n+k}$ ) 的 預測狀態。因此我們得到在時刻 $t_n, t_{n+1},...,t_{n+N}$ 對系統狀態的預測 (與狀態輸入序列 $u(0),...,u(N-1)$有關)。
每時刻的最佳控制力求解：
現在我們使用最佳控制的想法來決定 $u(0),...,u(N-1)$ 使得預測狀態 $x_u$ 能盡可能靠近 $x^{ref} = 0$。為了達成此目標，我們可選定 一距離函數
\[
l(x_u(k),u(k)) : = ||x_u(k)||^2 + \lambda ||u(k)||^2
\]其中 $||\cdot||$ 為 Euclidean norm 且 $\lambda \geq 0$為控制力的權重。則現在我們可以給出最佳控制問題如下
\[\min J\left( {x\left( n \right),u\left( . \right)} \right): = \sum\limits_{k = 0}^{N - 1} {l\left( {{x_u}\left( k \right),u\left( k \right)} \right)} \]對 $u(0),...,u(N-1)$ 為admissible 且 $x_u$ 滿足 (**)。
讓我們假設上述最佳控制問題具有解，且此解由一組最佳控制數列 記作
\[
u^*(0),...,u^*(N-1)
\]亦即我們有
\[\mathop {\min }\limits_{u\left( 0 \right),...,u\left( {N - 1} \right)} J\left( {x\left( n \right),u\left( . \right)} \right) = \sum\limits_{k = 0}^{N - 1} {l\left( {{x_{{u^*}}}\left( k \right),{u^*}\left( k \right)} \right)} \]現在為了讓控制力 $u$ 具有我們需要的回授形式 $\mu(x(n))$ 我們令
\[
\mu(x(n) ) := u^*(0)
\]亦即我們只使用最佳控制序列的第一元素。接著我們重複上述過程：在時刻 $t_{n+1}$ 可量得 $x(n+1)$ 並依此再度執行最佳化求得 $\mu(x(n+1))$。同理，對於時刻 $t_{n+2}$ 可量得 $x_{n+2}$ 並依此執行最佳化得到 $\mu(x(n+2))$...
Comments:
1. 因為受控系統為非線性，一般而言為我們所得到的 回授控制律 $\mu(\cdot)$ 必須透過迭代最佳化演算 而得，沒有解析解，故計算複雜度會是模型預測控制的一大挑戰。
2. 就預測模型觀點而言，我們可看出 預測狀態軌跡 $x_u(k), \;\; k=0,1,2,...,N$ 提供了對原本系統在已知時刻 $t_n$ 對於 $t_n,...,t_{n+N}$ 的預測，以及在時刻 $t_{n+1}$ 對於 $t_{n+1},...,t_{n+N+1}$ 的預測，以及 在時刻 $t_{n+2}$ 對於 $t_{n+2},...,t_{n+N+2}$ 的預測。因此我們可看出預測區間是 "移動" 的，故模型預測控制一般又稱為 移動區間控制 (Moving Horizon Control, MHC) 或稱 (Receding Horizon Control, RHC)
3. 若上述 $f$ 為線性，則上述討論稱之為 線性模型預測控制(Linear Model Predictive Control, LMPC)或者簡稱 模型預測控制 (Model Predictive Control, MPC)
4/06/2017
[訊號與系統] FIR系統的弦波響應 (1)
考慮 FIR 系統為 線性非時變(Linear Time-Invariant, LTI) 系統，且假設輸入為 離散 complex exponential 則其對應的輸出將非常容易計算：考慮 FIR 系統
\[
y[n] = \sum_{k=0}^M b_k x[n-k]
\]
假設輸入為 complex exponential 表為 $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 且 $-\infty < n < \infty$。則輸出為
\[\begin{array}{l}
y[n] = \sum\limits_{k = 0}^M {{b_k}} x[n - k]\\
= \sum\limits_{k = 0}^M {{b_k}} A{e^{j\varphi }}{e^{j\omega Ts\left( {n - k} \right)}}\\
= \left( {\sum\limits_{k = 0}^M {{b_k}} {e^{j\omega Ts\left( { - k} \right)}}} \right)A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}\\
:= H\left( {{e^{ - j\widehat \omega }}} \right)\underbrace {A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}}_{ = x\left[ n \right]}
\end{array}\]其中 $\widehat{\omega} := \omega T_s$ 且
\[H\left( {{e^{ - j\widehat \omega }}} \right) = \sum\limits_{k = 0}^M {{b_k}} {e^{ - j\widehat \omega k}} = \sum\limits_{k = 0}^M {h\left[ k \right]} {e^{ - j\widehat \omega k}}\]稱作 frequency-response function，一般而言我們簡稱為 frequency response。
回憶由於 FIR系統的脈衝響應 與 濾波器的係數相同，亦即 $b_k = h[k]$，我們可以將上述頻率響應改寫成
\[H\left( {{e^{ - j\widehat \omega k}}} \right) = \sum\limits_{k = 0}^M {{b_k}} {e^{ - j\widehat \omega k}} = \sum\limits_{k = 0}^M {h\left[ k \right]} {e^{ - j\widehat \omega k}}\]
Comments:
1. 輸入為離散時間訊號 $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 則 LTI FIR 系統亦為離散時間 complex exponential 乘上不同的複數大小，但頻率 $\widehat{\omega}$ 不變。
2. 上述以離散時間訊號為 complex exponential $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 為輸入，則 LTI FIR 系統可表為
\[y[n] = H\left( {{e^{ - j\widehat \omega k}}} \right)x\left[ n \right]\]
3. 注意到 $H\left( {{e^{ - j\widehat \omega k}}} \right)$ 為複數值函數，故我們可使用極座標表示法將其表示為
\[y[n] = \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|{e^{j\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}}} \right)A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}\]
亦即，
\[\begin{array}{l}
y[n] = H\left( {{e^{ - j\widehat \omega k}}} \right)\underbrace {A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}}_{ = x\left[ n \right]}\\
= \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|{e^{j\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}}} \right)A{e^{j\varphi }}{e^{j\widehat \omega n}}\\
= \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|A} \right){e^{j\left( {\angle H\left( {{e^{ - j\widehat \omega k}}} \right) + \varphi } \right)}}{e^{j\widehat \omega n}}
\end{array}\]一般而言，我們會稱 ${\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|}$ 為系統增益 (gain) 且 ${\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}$ 稱為系統相位 (phase)。
\[
y[n] = \sum_{k=0}^M b_k x[n-k]
\]
假設輸入為 complex exponential 表為 $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 且 $-\infty < n < \infty$。則輸出為
\[\begin{array}{l}
y[n] = \sum\limits_{k = 0}^M {{b_k}} x[n - k]\\
= \sum\limits_{k = 0}^M {{b_k}} A{e^{j\varphi }}{e^{j\omega Ts\left( {n - k} \right)}}\\
= \left( {\sum\limits_{k = 0}^M {{b_k}} {e^{j\omega Ts\left( { - k} \right)}}} \right)A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}\\
:= H\left( {{e^{ - j\widehat \omega }}} \right)\underbrace {A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}}_{ = x\left[ n \right]}
\end{array}\]其中 $\widehat{\omega} := \omega T_s$ 且
\[H\left( {{e^{ - j\widehat \omega }}} \right) = \sum\limits_{k = 0}^M {{b_k}} {e^{ - j\widehat \omega k}} = \sum\limits_{k = 0}^M {h\left[ k \right]} {e^{ - j\widehat \omega k}}\]稱作 frequency-response function，一般而言我們簡稱為 frequency response。
回憶由於 FIR系統的脈衝響應 與 濾波器的係數相同，亦即 $b_k = h[k]$，我們可以將上述頻率響應改寫成
\[H\left( {{e^{ - j\widehat \omega k}}} \right) = \sum\limits_{k = 0}^M {{b_k}} {e^{ - j\widehat \omega k}} = \sum\limits_{k = 0}^M {h\left[ k \right]} {e^{ - j\widehat \omega k}}\]
Comments:
1. 輸入為離散時間訊號 $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 則 LTI FIR 系統亦為離散時間 complex exponential 乘上不同的複數大小，但頻率 $\widehat{\omega}$ 不變。
2. 上述以離散時間訊號為 complex exponential $x[n] = x(nT_s) = A e^{j \varphi} e^{j \omega Ts n}$ 為輸入，則 LTI FIR 系統可表為
\[y[n] = H\left( {{e^{ - j\widehat \omega k}}} \right)x\left[ n \right]\]
3. 注意到 $H\left( {{e^{ - j\widehat \omega k}}} \right)$ 為複數值函數，故我們可使用極座標表示法將其表示為
\[y[n] = \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|{e^{j\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}}} \right)A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}\]
亦即，
\[\begin{array}{l}
y[n] = H\left( {{e^{ - j\widehat \omega k}}} \right)\underbrace {A{e^{j\varphi }}{e^{j\omega Ts\left( n \right)}}}_{ = x\left[ n \right]}\\
= \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|{e^{j\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}}} \right)A{e^{j\varphi }}{e^{j\widehat \omega n}}\\
= \left( {\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|A} \right){e^{j\left( {\angle H\left( {{e^{ - j\widehat \omega k}}} \right) + \varphi } \right)}}{e^{j\widehat \omega n}}
\end{array}\]一般而言，我們會稱 ${\left| {H\left( {{e^{ - j\widehat \omega k}}} \right)} \right|}$ 為系統增益 (gain) 且 ${\angle H\left( {{e^{ - j\widehat \omega k}}} \right)}$ 稱為系統相位 (phase)。
2/24/2017
[凸分析] 常見的凸集性質(1) - 線性矩陣不等式之解 所成的集合 為 凸集
給定 $a \in \mathbb{R}^n$ ，我們定義 線性函數 $f : \mathbb{R}^n \to \mathbb{R}$ 滿足
\[
f(x) := a^T x = a_1 x_1 + ... + a_n x_n
\]
現在我們進一步推廣上述結果：亦即上述的向量 $a = (a_1,...,a_n)$ 可以用 對稱矩陣 $(A_1,...,A_n)$ 替換，且 $ A_i \in S^m$ 為 $\mathbb{R}^{m \times m}$ 對稱矩陣，現在我們模仿上述線性函數 $f$ 定義一個新的函數如下：定義 $F: \mathbb{R}^n \to S^m$ 滿足
\[
F(x) := x_1 A_1 + ... + x_n A_n
\]
Comments:
1. $ F(x) $ 仍為 $\mathbb{R}^{m \times m}$ 的對稱矩陣。
2. 上述提及的 線性函數 $f(x)$ (或者又說標準內積 或者 hyperplane) 可用以形成所謂 convex polyhedra 的集合，在此不贅述。
接著我們想問 對於上述 矩陣等式 $g(x)$ 而言，是否可以定義不等式？ 一般而言在線性代數中我們定義 $F(x) \succ 0$ 表示 $F(x)$ 為正定矩陣，亦即 對任意 $z \in \mathbb{R}^n$ 且 $z \neq 0$ 我們有
\[
z^T F(x) z > 0
\] 我們說 $F(x) \succeq 0$ 表示 $F(x)$ 為半正定矩陣，亦即 對任意 $z \in \mathbb{R}^n$
\[
z^T F(x) z \geq 0
\]
FACT:
令 $A,B$ 為 兩實係數 對稱矩陣，若 $A \succeq 0$ 且 $B \succeq 0$ 則
\[
A+B \succeq 0
\]
Proof: omitted (此證明相對容易，在此略過)
========================
Definition: Linear Matrix Inequality (LMI)
我們稱一不等式 為對 $x$ 而言的線性矩陣不等式 (Linear Matrix Inequality in $x$, LMI) 若 前述的矩陣 $F(x)$ 具有下列形式：
\[
F(x) := x_1 A_1 + ... + x_n A_n \preceq B
\] 其中 $x_i \in \mathbb{R}^1$ 且 $B, A_i $ 為 $m \times m$ 對稱矩陣，$i=1,2,...,n$。
========================
Comments:
1. 上述 LMI 要求 $F(x) \preceq B $ 亦即 $B - F(x) \succeq 0$ ，也就是說 $B - F(x) $ 為正定對稱矩陣，由前述定義可知我們要求：對任意 $z \in \mathbb{R}^n$，
\[
z^T (B-F(x))z \geq 0
\]
2. LMI 為 "線性" in $x$
3. LMI 在 強健控制理論中扮演重要的角色，在此不贅述。
以下我們給出主要結果：
========================
FACT:
上述 LMI 之解所成之集合 \[
L:=\{x \in \mathbb{R}^n : F(x) \preceq B\}
\]為凸集。
========================
Proof:
令 $x,y \in L$ 且 $\theta \in [0,1]$ 我們要證明 $ \theta x + (1-\theta)y \in L $ 此等價於證明
\[
F(\theta x + (1-\theta)y) \preceq B
\] 現在觀察
\begin{align*}
F(\theta x + (1 - \theta )y) &= (\theta {x_1} + (1 - \theta ){y_1}){A_1} + ... + (\theta {x_n} + (1 - \theta ){y_n}){A_n} \hfill \\
&= \sum\limits_{i = 1}^n {(\theta {x_i} + (1 - \theta ){y_i}){A_i}} \hfill \\
&= \theta \sum\limits_{i = 1}^n {{x_i}{A_i}} + (1 - \theta )\sum\limits_{i = 1}^n {{y_i}{A_i}} \;\;\;\;\; (*) \hfill \\
\end{align*}
由於 $x,y \in L$ ，故我們有
\begin{align*}
F(x) &:= \sum_{i=1}^n x_i A_i \preceq B; \\
F(y) &:= \sum_{i=1}^n y_i A_i \preceq B
\end{align*}故將上述結果帶入 $(*)$ ，由於 $\theta \in [0,1]$ 利用前述 FACT 可得
\[
F(\theta x + (1-\theta)y) \preceq B
\]至此得證。$\square$
\[
f(x) := a^T x = a_1 x_1 + ... + a_n x_n
\]
現在我們進一步推廣上述結果：亦即上述的向量 $a = (a_1,...,a_n)$ 可以用 對稱矩陣 $(A_1,...,A_n)$ 替換，且 $ A_i \in S^m$ 為 $\mathbb{R}^{m \times m}$ 對稱矩陣，現在我們模仿上述線性函數 $f$ 定義一個新的函數如下：定義 $F: \mathbb{R}^n \to S^m$ 滿足
\[
F(x) := x_1 A_1 + ... + x_n A_n
\]
Comments:
1. $ F(x) $ 仍為 $\mathbb{R}^{m \times m}$ 的對稱矩陣。
2. 上述提及的 線性函數 $f(x)$ (或者又說標準內積 或者 hyperplane) 可用以形成所謂 convex polyhedra 的集合，在此不贅述。
接著我們想問 對於上述 矩陣等式 $g(x)$ 而言，是否可以定義不等式？ 一般而言在線性代數中我們定義 $F(x) \succ 0$ 表示 $F(x)$ 為正定矩陣，亦即 對任意 $z \in \mathbb{R}^n$ 且 $z \neq 0$ 我們有
\[
z^T F(x) z > 0
\] 我們說 $F(x) \succeq 0$ 表示 $F(x)$ 為半正定矩陣，亦即 對任意 $z \in \mathbb{R}^n$
\[
z^T F(x) z \geq 0
\]
FACT:
令 $A,B$ 為 兩實係數 對稱矩陣，若 $A \succeq 0$ 且 $B \succeq 0$ 則
\[
A+B \succeq 0
\]
Proof: omitted (此證明相對容易，在此略過)
========================
Definition: Linear Matrix Inequality (LMI)
我們稱一不等式 為對 $x$ 而言的線性矩陣不等式 (Linear Matrix Inequality in $x$, LMI) 若 前述的矩陣 $F(x)$ 具有下列形式：
\[
F(x) := x_1 A_1 + ... + x_n A_n \preceq B
\] 其中 $x_i \in \mathbb{R}^1$ 且 $B, A_i $ 為 $m \times m$ 對稱矩陣，$i=1,2,...,n$。
========================
1. 上述 LMI 要求 $F(x) \preceq B $ 亦即 $B - F(x) \succeq 0$ ，也就是說 $B - F(x) $ 為正定對稱矩陣，由前述定義可知我們要求：對任意 $z \in \mathbb{R}^n$，
\[
z^T (B-F(x))z \geq 0
\]
2. LMI 為 "線性" in $x$
3. LMI 在 強健控制理論中扮演重要的角色，在此不贅述。
以下我們給出主要結果：
========================
FACT:
上述 LMI 之解所成之集合 \[
L:=\{x \in \mathbb{R}^n : F(x) \preceq B\}
\]為凸集。
========================
令 $x,y \in L$ 且 $\theta \in [0,1]$ 我們要證明 $ \theta x + (1-\theta)y \in L $ 此等價於證明
\[
F(\theta x + (1-\theta)y) \preceq B
\] 現在觀察
\begin{align*}
F(\theta x + (1 - \theta )y) &= (\theta {x_1} + (1 - \theta ){y_1}){A_1} + ... + (\theta {x_n} + (1 - \theta ){y_n}){A_n} \hfill \\
&= \sum\limits_{i = 1}^n {(\theta {x_i} + (1 - \theta ){y_i}){A_i}} \hfill \\
&= \theta \sum\limits_{i = 1}^n {{x_i}{A_i}} + (1 - \theta )\sum\limits_{i = 1}^n {{y_i}{A_i}} \;\;\;\;\; (*) \hfill \\
\end{align*}
由於 $x,y \in L$ ，故我們有
\begin{align*}
F(x) &:= \sum_{i=1}^n x_i A_i \preceq B; \\
F(y) &:= \sum_{i=1}^n y_i A_i \preceq B
\end{align*}故將上述結果帶入 $(*)$ ，由於 $\theta \in [0,1]$ 利用前述 FACT 可得
\[
F(\theta x + (1-\theta)y) \preceq B
\]至此得證。$\square$
8/11/2016
[變分法] 連續泛函極值的必要條件
這次要介紹最簡單形式的 泛函極值問題的 必要條件，此條件一般又稱之為 Euler-Largrange Eqution。此方程可謂泛函極值的房角石，亦為之後在最佳控制理論中的最大值原理扮演開路先鋒，是極為重要的角色。在介紹之前，我們先做一般性的用語與基本性質介紹。
======================
Definition: 泛函
令 $\Omega$ 為 賦範函數空間 (normed function space)，若 對任意函數 $x(t) \in \Omega$ 都存在一個實數與之對應，則我們稱 $J$ 是定義在 $\Omega$ 上的 泛函 (functional)，記作 $J(x(t))$
======================
Comment:
1. 簡而言之，泛函 一詞即表示為由 函數空間 映射到 實數軸 上的函數 $J: \Omega \to \mathbb{R}$ 。
2. 再以下的討論中，集合 $\Omega$ 又稱為 泛函 $J$ 的 容許集 (admissible set)。
現取 $x_1, x \in \Omega$ 且 $\delta x := x_1 - x$ ，我們定義 關於 $\delta x$ 的 泛函增量 (increment) 如下
\[
\Delta J(\delta x) := J(x_1) - J(x) = J( x + \delta x) - J(x)
\]則由此 泛函增量，我們可以定義何謂泛函的變分。
======================
Definition: 泛函的變分
給定泛函 $J : \Omega \to \mathbb{R}$，若存在 一線性泛函 $L(x, \delta x)$ 使得泛函增量可被表為
\[
\Delta J(\delta x) = L(x, \delta x) + r(x, \delta x) \cdot | |\delta x||
\]其中 $r(x, \delta x)$ 為 其他高階剩餘項(remainder) 滿足 當 $| |\delta x|| \to 0 \Rightarrow r(x, \delta x) \to 0$，則我們稱上式中的 $L(x, \delta x)$ 為 $J(x)$ 的 變分 (variation)，記作 $\delta J := L(x, \delta x)$
======================
======================
Theorem: 泛函極值與變分關係
給定泛函 $J : \Omega \to \mathbb{R}$，若其變分存在，則 其變分可表為參數 $\alpha$ 的方向導數，亦即 變分滿足下式
\[\delta J(x(t)) = {\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}}\]======================
Proof: 給定泛函 $J : \Omega \to \mathbb{R}$ 且假設其變分存在，我們要證明
\[\delta J(x(t)) = {\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}}
\]首先由 $\delta J$ 存在可知：存在一線性泛函 $L$ 始得 泛函增量 $\Delta J$滿足
\[\begin{align*}
\Delta J &= J\left( {x + \alpha \delta x} \right) - J\left( x \right) \hfill \\
&= L(x,\alpha \delta x) + r(x,\alpha \delta x) \cdot || \alpha \delta x || \hfill \\
\end{align*}
\]由於 $L$ 為線性泛函，故 $L(x,\alpha \delta x) = \alpha L(x,\delta x)$，現在觀察
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}} &= \mathop {\lim }\limits_{\alpha \to 0} \frac{{J(x + \alpha \delta x) - J\left( x \right)}}{\alpha }\\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{L(x,\alpha \delta x) + r(x,\alpha \delta x)||\alpha \delta x||}}{\alpha } \hfill \\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{L(x,\alpha \delta x)}}{\alpha } + \underbrace {\mathop {\lim }\limits_{\alpha \to 0} \frac{{r(x,\alpha \delta x) ||\alpha \delta x||}}{\alpha }}_{ = 0} \hfill \\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{\alpha L(x,\delta x)}}{\alpha } \hfill \\
&= L(x,\delta x) = \delta J(x) \;\;\;\;\; \square
\end{align*}
\]
======================
Theorem:
令 $J$ 為泛函且其變分存在，若 $J(x)$ 在 $x_0 \in \Omega$ 有(局部)極值，則其在 $x_0$ 之變分
\[
\delta J(x_0) =0
\] ======================
\[{\left. {\delta J\left( x \right) = \frac{\partial }{{\partial \alpha }}J(x + \alpha \delta x)} \right|_{\alpha = 0}}\]由於 $J(x)$ 在 $x_0 \in \Omega$ 有局部極值，故我們可知 $\alpha =0$ 為$J(x_0 + \alpha \delta x)$ 的局部極值 (以極小值為例，可知對任意 $\alpha \in \mathbb{R}$， $J(x_0) \leq J(x_0 + \alpha \delta x)$，且極小值發生在 $\alpha = 0$)，故
\[{\left. {\delta J\left( {{x_0}} \right) = \frac{\partial }{{\partial \alpha }}J({x_0} + \alpha \delta x)} \right|_{\alpha = 0}} = 0\;\;\;\; \square
\]
在討論一般設定之後，以下我們開始針對特殊形式的泛函來建構必要條件：考慮泛函
\[
J(x(t)) := \int_{t_0}^{t_1} F(t,x,\dot{x}) dt; \;\;\; x(t_0) :=x_0; \;\;\; x(t_1) \doteq x_1
\]且令其 admissible set 為
\[
\Omega := \{x(t) : x(t) \in C^2[t_0, t_1], \; x(t_0) = x_0, x(t_1) = x_1\}
\]且 $F(t, x, \dot{x})$ 為 $C^2$ (二階可導且連續)，我們欲求上述泛函極值的必要條件，此結果極為鼎鼎大名的 Euler-Largrange 方程，但在我們證明主要定理之前，底下我們先給個前置定理，此定理又稱為變分基本定理。
======================
Lemma: 變分基本引理
設函數 $F(t)$ 在區間 $[t_0, t_1]$ 上連續，若對於任意滿足 $\eta(t_0) = \eta(t_1) =0$ 的充分光滑函數 $\eta(t)$ 我們都有
\[
\int_{t_0}^{t_1} F(t) \eta(t) dt =0
\]則 $F(t) = 0$ 對 $t \in [t_0,t_1]$
======================
Proof: 利用反證法，假設 存在 $\xi \in (t_0,t_1)$ 使得 $F(\xi) \neq 0$，欲證明矛盾。不失一般性情況下我們假設 $F(\xi) >0$ 則由於 $F$的連續性，可知必存在 以 $\xi$ 為中心的鄰域 $N_\xi :=(\xi_1,\xi_2) \subset (t_0, t_1)$ 使得 對任意 $t \in N_\xi$，我們有 $F(t) > 0$。現在我們構造 $\eta(t)$ 函數如下
\[\eta \left( t \right): = \left\{ \begin{gathered}
0,\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}t \in \left[ {{t_0},{\xi _1}} \right) \hfill \\
{\left[ {\left( {t - {\xi _1}} \right)\left( {t - {\xi _2}} \right)} \right]^2},\begin{array}{*{20}{c}}
{}&{}
\end{array}t \in \left[ {{\xi _1},{\xi _2}} \right] \hfill \\
0,\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}t \in \left( {{\xi _2},{t_1}} \right] \hfill \\
\end{gathered} \right.\]且注意到上述 $\eta(t)$ 函數滿足 $\eta(t_0) = \eta(t_1) = 0$ 且為連續函數，然而若我們觀察
\[
\int_{t_0}^{t_1} F(t) \eta(t) dt = \int_{\xi_1}^{\xi_2} F(t) \eta(t) dt > 0
\]此結果與我們的假設矛盾。$\square$
======================
Theorem: 泛函極值的必要條件 Euler-Lagrange Equation
設函數 $F(t, x, \dot{x})$ 具有連續二階偏導數，且設泛函\[
J(x(t)) := \int_{t_0}^{t_1} F(t,x,\dot{x}) dt; \;\;\; x(t_0) :=x_0; \;\;\; x(t_1) \doteq x_1
\]在 $x(t) \in \Omega$ 達到極值，則 $x(t)$ 滿足下列方程
\[\frac{\partial }{{\partial x}}F\left( {t,x,\dot x} \right) - \frac{d}{{dt}}\left( {\frac{\partial }{{\partial \dot x}}F\left( {t,x,\dot x} \right)} \right) = 0\]
======================
Proof: 首先令 $\phi(t) := \delta x(t)$，則由於 $x(t_0)=x_0$與 $x(t_1) = x_1$ 可知，$\phi(t)$ 滿足 $\phi(t_0) = \phi(t_1)=0$，現在由泛函極值與變分的關係可知下式必定成立：
\[\delta J\left( {x\left( t \right)} \right) = \left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0} = 0 \;\;\;(\star)
\]現在觀察
\[
J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right) = \int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt
\]故我們可先行計算
\[{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}} = {\left. {\frac{\partial }{{\partial \alpha }}\int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}}
\]由 Libneiz Rule 可得
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}} &= {\left. {\frac{\partial }{{\partial \alpha }}\int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\frac{\partial }{{\partial \alpha }}F} (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi + \frac{{\partial F}}{{\partial \dot x}}\dot \phi } \right]} dt} \right|_{\alpha = 0}}\;\;\;\; (*)
\end{align*}
\]注意到上述積分第二項可透過 integration by part 求得
\[\int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = \left. {\frac{{\partial F}}{{\partial \dot x}}\phi } \right|_{{t_0}}^{{t_1}} - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt\]由於 $\phi(t_0) = \phi(t_1) = 0$，故我們得
\[\begin{gathered}
\int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = \underbrace {\left. {\frac{{\partial F}}{{\partial \dot x}}\phi } \right|_{{t_0}}^{{t_1}}}_{ = 0} - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt \hfill \\
\Rightarrow \int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt \hfill \\
\end{gathered}
\]現在將其帶回 $(*)$ 我們得到
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}}
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi - \phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} \phi dt} \right|_{\alpha = 0}} \hfill \\
\end{align*}
\]由於 $(\star)$ 可知，
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}}
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi + \frac{{\partial F}}{{\partial \dot x}}\dot \phi } \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi - \phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} \phi dt} \right|_{\alpha = 0}} = 0 \hfill \\
\end{align*}
\]由於 ${\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}}$ 在區間 $[t_0,t_1]$ 連續，且 $\phi$ 滿足 $\phi(t_0) = \phi(t_1) =0$ 且 $\phi \in C^2$，利用前述引理可知在 $[t_0,t_1]$ 上，
\[{\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} = 0\;\;\;\; \square\]
[1] I. M. Gelfand and S. V. Fomin, Calculus of Variations, 2000
[2] David G. Luenberger, Optimization By Vector Space Methods, 1997
======================
Definition: 泛函
令 $\Omega$ 為 賦範函數空間 (normed function space)，若 對任意函數 $x(t) \in \Omega$ 都存在一個實數與之對應，則我們稱 $J$ 是定義在 $\Omega$ 上的 泛函 (functional)，記作 $J(x(t))$
======================
1. 簡而言之，泛函 一詞即表示為由 函數空間 映射到 實數軸 上的函數 $J: \Omega \to \mathbb{R}$ 。
2. 再以下的討論中，集合 $\Omega$ 又稱為 泛函 $J$ 的 容許集 (admissible set)。
現取 $x_1, x \in \Omega$ 且 $\delta x := x_1 - x$ ，我們定義 關於 $\delta x$ 的 泛函增量 (increment) 如下
\[
\Delta J(\delta x) := J(x_1) - J(x) = J( x + \delta x) - J(x)
\]則由此 泛函增量，我們可以定義何謂泛函的變分。
======================
Definition: 泛函的變分
給定泛函 $J : \Omega \to \mathbb{R}$，若存在 一線性泛函 $L(x, \delta x)$ 使得泛函增量可被表為
\[
\Delta J(\delta x) = L(x, \delta x) + r(x, \delta x) \cdot | |\delta x||
\]其中 $r(x, \delta x)$ 為 其他高階剩餘項(remainder) 滿足 當 $| |\delta x|| \to 0 \Rightarrow r(x, \delta x) \to 0$，則我們稱上式中的 $L(x, \delta x)$ 為 $J(x)$ 的 變分 (variation)，記作 $\delta J := L(x, \delta x)$
======================
Comment:
1. 上述定義中的 線性泛函項 $L$ 與 其他高階剩餘項 $r$，可視為透過 Taylor 級數展開而得。
2. 變分 (variation) 一詞在文獻中又稱 differential
3. 若泛函變分存在，則該 變分 為唯一，在此不證明，有興趣讀者可參閱 [1]。
4. 關於線性泛函及其相關定義請讀者可參閱 [變分法] 淺論 線性泛函
5. 有些文獻定義的泛函是透過所謂 Gateaux differentials 與 Freshet differential，但為求論述簡潔，在此不多作介紹，有興趣的讀者可以參閱 [2]
4. 關於線性泛函及其相關定義請讀者可參閱 [變分法] 淺論 線性泛函
5. 有些文獻定義的泛函是透過所謂 Gateaux differentials 與 Freshet differential，但為求論述簡潔，在此不多作介紹，有興趣的讀者可以參閱 [2]
======================
Theorem: 泛函極值與變分關係
給定泛函 $J : \Omega \to \mathbb{R}$，若其變分存在，則 其變分可表為參數 $\alpha$ 的方向導數，亦即 變分滿足下式
\[\delta J(x(t)) = {\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}}\]======================
\[\delta J(x(t)) = {\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}}
\]首先由 $\delta J$ 存在可知：存在一線性泛函 $L$ 始得 泛函增量 $\Delta J$滿足
\[\begin{align*}
\Delta J &= J\left( {x + \alpha \delta x} \right) - J\left( x \right) \hfill \\
&= L(x,\alpha \delta x) + r(x,\alpha \delta x) \cdot || \alpha \delta x || \hfill \\
\end{align*}
\]由於 $L$ 為線性泛函，故 $L(x,\alpha \delta x) = \alpha L(x,\delta x)$，現在觀察
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J(x(t) + \alpha \delta x(t))} \right|_{\alpha = 0}} &= \mathop {\lim }\limits_{\alpha \to 0} \frac{{J(x + \alpha \delta x) - J\left( x \right)}}{\alpha }\\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{L(x,\alpha \delta x) + r(x,\alpha \delta x)||\alpha \delta x||}}{\alpha } \hfill \\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{L(x,\alpha \delta x)}}{\alpha } + \underbrace {\mathop {\lim }\limits_{\alpha \to 0} \frac{{r(x,\alpha \delta x) ||\alpha \delta x||}}{\alpha }}_{ = 0} \hfill \\
&= \mathop {\lim }\limits_{\alpha \to 0} \frac{{\alpha L(x,\delta x)}}{\alpha } \hfill \\
&= L(x,\delta x) = \delta J(x) \;\;\;\;\; \square
\end{align*}
\]
======================
Theorem:
令 $J$ 為泛函且其變分存在，若 $J(x)$ 在 $x_0 \in \Omega$ 有(局部)極值，則其在 $x_0$ 之變分
\[
\delta J(x_0) =0
\] ======================
Comment: 上述定理中的 $x_0$ 又稱為 泛函 $J$ 的臨界點(critical point) 或者稱 不動點 (stationary point)。
Proof: 由於變分存在，我們可將變分用 以單變數參數 $\alpha$ 的方向導數表示\[{\left. {\delta J\left( x \right) = \frac{\partial }{{\partial \alpha }}J(x + \alpha \delta x)} \right|_{\alpha = 0}}\]由於 $J(x)$ 在 $x_0 \in \Omega$ 有局部極值，故我們可知 $\alpha =0$ 為$J(x_0 + \alpha \delta x)$ 的局部極值 (以極小值為例，可知對任意 $\alpha \in \mathbb{R}$， $J(x_0) \leq J(x_0 + \alpha \delta x)$，且極小值發生在 $\alpha = 0$)，故
\[{\left. {\delta J\left( {{x_0}} \right) = \frac{\partial }{{\partial \alpha }}J({x_0} + \alpha \delta x)} \right|_{\alpha = 0}} = 0\;\;\;\; \square
\]
在討論一般設定之後，以下我們開始針對特殊形式的泛函來建構必要條件：考慮泛函
\[
J(x(t)) := \int_{t_0}^{t_1} F(t,x,\dot{x}) dt; \;\;\; x(t_0) :=x_0; \;\;\; x(t_1) \doteq x_1
\]且令其 admissible set 為
\[
\Omega := \{x(t) : x(t) \in C^2[t_0, t_1], \; x(t_0) = x_0, x(t_1) = x_1\}
\]且 $F(t, x, \dot{x})$ 為 $C^2$ (二階可導且連續)，我們欲求上述泛函極值的必要條件，此結果極為鼎鼎大名的 Euler-Largrange 方程，但在我們證明主要定理之前，底下我們先給個前置定理，此定理又稱為變分基本定理。
======================
Lemma: 變分基本引理
設函數 $F(t)$ 在區間 $[t_0, t_1]$ 上連續，若對於任意滿足 $\eta(t_0) = \eta(t_1) =0$ 的充分光滑函數 $\eta(t)$ 我們都有
\[
\int_{t_0}^{t_1} F(t) \eta(t) dt =0
\]則 $F(t) = 0$ 對 $t \in [t_0,t_1]$
======================
\[\eta \left( t \right): = \left\{ \begin{gathered}
0,\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}t \in \left[ {{t_0},{\xi _1}} \right) \hfill \\
{\left[ {\left( {t - {\xi _1}} \right)\left( {t - {\xi _2}} \right)} \right]^2},\begin{array}{*{20}{c}}
{}&{}
\end{array}t \in \left[ {{\xi _1},{\xi _2}} \right] \hfill \\
0,\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}t \in \left( {{\xi _2},{t_1}} \right] \hfill \\
\end{gathered} \right.\]且注意到上述 $\eta(t)$ 函數滿足 $\eta(t_0) = \eta(t_1) = 0$ 且為連續函數，然而若我們觀察
\[
\int_{t_0}^{t_1} F(t) \eta(t) dt = \int_{\xi_1}^{\xi_2} F(t) \eta(t) dt > 0
\]此結果與我們的假設矛盾。$\square$
======================
Theorem: 泛函極值的必要條件 Euler-Lagrange Equation
設函數 $F(t, x, \dot{x})$ 具有連續二階偏導數，且設泛函\[
J(x(t)) := \int_{t_0}^{t_1} F(t,x,\dot{x}) dt; \;\;\; x(t_0) :=x_0; \;\;\; x(t_1) \doteq x_1
\]在 $x(t) \in \Omega$ 達到極值，則 $x(t)$ 滿足下列方程
\[\frac{\partial }{{\partial x}}F\left( {t,x,\dot x} \right) - \frac{d}{{dt}}\left( {\frac{\partial }{{\partial \dot x}}F\left( {t,x,\dot x} \right)} \right) = 0\]
======================
\[\delta J\left( {x\left( t \right)} \right) = \left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0} = 0 \;\;\;(\star)
\]現在觀察
\[
J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right) = \int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt
\]故我們可先行計算
\[{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}} = {\left. {\frac{\partial }{{\partial \alpha }}\int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}}
\]由 Libneiz Rule 可得
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}} &= {\left. {\frac{\partial }{{\partial \alpha }}\int_{{t_0}}^{{t_1}} F (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\frac{\partial }{{\partial \alpha }}F} (t,x + \alpha \phi ,\dot x + \alpha \dot \phi )dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi + \frac{{\partial F}}{{\partial \dot x}}\dot \phi } \right]} dt} \right|_{\alpha = 0}}\;\;\;\; (*)
\end{align*}
\]注意到上述積分第二項可透過 integration by part 求得
\[\int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = \left. {\frac{{\partial F}}{{\partial \dot x}}\phi } \right|_{{t_0}}^{{t_1}} - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt\]由於 $\phi(t_0) = \phi(t_1) = 0$，故我們得
\[\begin{gathered}
\int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = \underbrace {\left. {\frac{{\partial F}}{{\partial \dot x}}\phi } \right|_{{t_0}}^{{t_1}}}_{ = 0} - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt \hfill \\
\Rightarrow \int_{{t_0}}^{{t_1}} {\frac{{\partial F}}{{\partial \dot x}}\dot \phi dt} = - \int_{{t_0}}^{{t_1}} {\phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} dt \hfill \\
\end{gathered}
\]現在將其帶回 $(*)$ 我們得到
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}}
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi - \phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} \phi dt} \right|_{\alpha = 0}} \hfill \\
\end{align*}
\]由於 $(\star)$ 可知，
\[\begin{align*}
{\left. {\frac{\partial }{{\partial \alpha }}J\left( {x\left( t \right) + \alpha \phi \left( t \right)} \right)} \right|_{\alpha = 0}}
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi + \frac{{\partial F}}{{\partial \dot x}}\dot \phi } \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}}\phi - \phi \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} dt} \right|_{\alpha = 0}} \hfill \\
&= {\left. {\int_{{t_0}}^{{t_1}} {\left[ {\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} \right]} \phi dt} \right|_{\alpha = 0}} = 0 \hfill \\
\end{align*}
\]由於 ${\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}}$ 在區間 $[t_0,t_1]$ 連續，且 $\phi$ 滿足 $\phi(t_0) = \phi(t_1) =0$ 且 $\phi \in C^2$，利用前述引理可知在 $[t_0,t_1]$ 上，
\[{\frac{{\partial F}}{{\partial x}} - \frac{d}{{dt}}\frac{{\partial F}}{{\partial \dot x}}} = 0\;\;\;\; \square\]
[1] I. M. Gelfand and S. V. Fomin, Calculus of Variations, 2000
[2] David G. Luenberger, Optimization By Vector Space Methods, 1997
3/28/2016
[控制理論] 具有負實部特徵值之 LTV系統 並不保證系統穩定
考慮線性非時變 (Linear Time-Invariant, LTI)系統利用狀態空間表示：
\[
{\bf \dot x} = A {\bf x} + B {\bf u}
\] 回憶在大學部自動控制課程中，我們知道 LTI 系統穩定 的 充分必要條件 為系統矩陣 $A$ 之特徵值具有負實部 （或者等價論述為 極點 pole 落在 複數平面的左半面）。現在我們想問若 系統為 線性時變 (Linear Time-Varying, LTV)系統是否此條件依然成立？
答案是否定的，以下為一個極為出色的反例：考慮線性時變系統 ${\bf \dot x} = A(t) {\bf x} $ 其中
\[A\left( t \right): = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]
\] 且給定初始狀態為 $x_1(0) = x_2(0)=1$ 則由於此系統 $A(t)$ 矩陣為三角矩陣，其特徵值為對角線元素，亦即 $\lambda_{1,2} = -1$，具有負實部。然而，若我們求解此 LTV 系統，亦即觀察
\[{\bf{\dot x}}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]\left[ \begin{array}{l}
{x_1}\left( t \right)\\
{x_2}\left( t \right)
\end{array} \right] = \left[ \begin{array}{l}
- {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)\\
- {x_2}\left( t \right)
\end{array} \right]\]故我們可首先解得
\[\begin{array}{*{20}{l}}
{{{\dot x}_2}\left( t \right) = - {x_2}\left( t \right)}\\
\begin{array}{l}
\Rightarrow {x_2}\left( t \right) = {e^{ - t}}{x_2}\left( 0 \right)\\
\Rightarrow {x_2}\left( t \right) = {e^{ - t}}
\end{array}
\end{array}
\]再將此 $x_2(t)$ 帶回 $\dot x_1(t)$ 式中，可求解 $x_1$ 如下
\[\begin{array}{*{20}{l}}
{{{\dot x}_1}\left( t \right) = - {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)}\\
{ \Rightarrow {{\dot x}_1}\left( t \right) = - {x_1}\left( t \right) + {e^{2t}}{e^{ - t}}}\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}}{x_1}\left( 0 \right) + \int_0^t {{e^{ - \left( {t - \tau } \right)}}{e^\tau }d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}} + {e^{ - \left( t \right)}}\int_0^t {{e^{2\tau }}d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = \frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}\]故系統之解為
\[{{\bf{x}}\left( t \right) = \left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]}\]但注意到若我們計算上述之狀態的 2-norm 且取極限 $t \to \infty$ 會發現
\[\begin{array}{l}
\mathop {\lim }\limits_{t \to \infty } \left\| {{\bf{x}}\left( t \right)} \right\| = \mathop {\lim }\limits_{t \to \infty } \left\| {\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right\|\\
= \mathop {\lim }\limits_{t \to \infty } {\left( {\left[ {\begin{array}{*{20}{c}}
{{e^{ - t}}}&{\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}} \right]\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right)^{1/2}}\\
= \mathop {\lim }\limits_{t \to \infty } {\left( {{e^{ - 2t}} + {{\left( {\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}} \right)}^2}} \right)^{1/2}} = \infty
\end{array}\]亦即系統狀態發散。
上述結果闡釋了對於 LTV 系統而言，負實部特徵值 (左半面極點) 不保證系統穩定。
\[
{\bf \dot x} = A {\bf x} + B {\bf u}
\] 回憶在大學部自動控制課程中，我們知道 LTI 系統穩定 的 充分必要條件 為系統矩陣 $A$ 之特徵值具有負實部 （或者等價論述為 極點 pole 落在 複數平面的左半面）。現在我們想問若 系統為 線性時變 (Linear Time-Varying, LTV)系統是否此條件依然成立？
答案是否定的，以下為一個極為出色的反例：考慮線性時變系統 ${\bf \dot x} = A(t) {\bf x} $ 其中
\[A\left( t \right): = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]
\] 且給定初始狀態為 $x_1(0) = x_2(0)=1$ 則由於此系統 $A(t)$ 矩陣為三角矩陣，其特徵值為對角線元素，亦即 $\lambda_{1,2} = -1$，具有負實部。然而，若我們求解此 LTV 系統，亦即觀察
\[{\bf{\dot x}}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{ - 1}&{{e^{2t}}}\\
0&{ - 1}
\end{array}} \right]\left[ \begin{array}{l}
{x_1}\left( t \right)\\
{x_2}\left( t \right)
\end{array} \right] = \left[ \begin{array}{l}
- {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)\\
- {x_2}\left( t \right)
\end{array} \right]\]故我們可首先解得
\[\begin{array}{*{20}{l}}
{{{\dot x}_2}\left( t \right) = - {x_2}\left( t \right)}\\
\begin{array}{l}
\Rightarrow {x_2}\left( t \right) = {e^{ - t}}{x_2}\left( 0 \right)\\
\Rightarrow {x_2}\left( t \right) = {e^{ - t}}
\end{array}
\end{array}
\]再將此 $x_2(t)$ 帶回 $\dot x_1(t)$ 式中，可求解 $x_1$ 如下
\[\begin{array}{*{20}{l}}
{{{\dot x}_1}\left( t \right) = - {x_1}\left( t \right) + {e^{2t}}{x_2}\left( t \right)}\\
{ \Rightarrow {{\dot x}_1}\left( t \right) = - {x_1}\left( t \right) + {e^{2t}}{e^{ - t}}}\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}}{x_1}\left( 0 \right) + \int_0^t {{e^{ - \left( {t - \tau } \right)}}{e^\tau }d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = {e^{ - t}} + {e^{ - \left( t \right)}}\int_0^t {{e^{2\tau }}d\tau } }\\
{ \Rightarrow {x_1}\left( t \right) = \frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}\]故系統之解為
\[{{\bf{x}}\left( t \right) = \left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]}\]但注意到若我們計算上述之狀態的 2-norm 且取極限 $t \to \infty$ 會發現
\[\begin{array}{l}
\mathop {\lim }\limits_{t \to \infty } \left\| {{\bf{x}}\left( t \right)} \right\| = \mathop {\lim }\limits_{t \to \infty } \left\| {\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right\|\\
= \mathop {\lim }\limits_{t \to \infty } {\left( {\left[ {\begin{array}{*{20}{c}}
{{e^{ - t}}}&{\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}}
\end{array}} \right]\left[ \begin{array}{l}
{e^{ - t}}\\
\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}
\end{array} \right]} \right)^{1/2}}\\
= \mathop {\lim }\limits_{t \to \infty } {\left( {{e^{ - 2t}} + {{\left( {\frac{1}{2}{e^t} + \frac{1}{2}{e^{ - \left( t \right)}}} \right)}^2}} \right)^{1/2}} = \infty
\end{array}\]亦即系統狀態發散。
上述結果闡釋了對於 LTV 系統而言，負實部特徵值 (左半面極點) 不保證系統穩定。
2/16/2016
[投資理論] 數學能擊敗金融市場嗎？
以下為個人在 University of Wisconsin-Madison 臺灣學生會 2016年 第一場學術沙龍中 分享的簡報：
數學能擊敗金融市場嗎？-從控制理論觀點 （2016, 02. 16)
講者：謝宗翰
講題：數學是否能擊敗金融市場？-從控制理論觀點
簡介：此講題將試圖回答一個基本問題：是否存在一種「必勝法」，使得投資績效具備恆正報酬？我們將從現代投資理論出發，最終止於財務工程與控制理論，過程中，我們將逐步揭示何時可以透過數學幫助我們建構一組可行的「最佳」交易策略。
關於 UW-Madison 臺灣學生會 連結
https://sites.google.com/site/satuwmadison/
數學能擊敗金融市場嗎？-從控制理論觀點 （2016, 02. 16)
講者：謝宗翰
講題：數學是否能擊敗金融市場？-從控制理論觀點
簡介：此講題將試圖回答一個基本問題：是否存在一種「必勝法」，使得投資績效具備恆正報酬？我們將從現代投資理論出發，最終止於財務工程與控制理論，過程中，我們將逐步揭示何時可以透過數學幫助我們建構一組可行的「最佳」交易策略。
關於 UW-Madison 臺灣學生會 連結
https://sites.google.com/site/satuwmadison/
9/27/2015
[自動控制] 穩態誤差與特性方程反求 轉移函數問題
考慮單位回授控制系統如下圖
現在假設
1. 閉迴路系統對 單位步階訊號 的穩態誤差為零：
2. 閉迴路轉移函數 $Y(s)/R(s)$ 的特性方程式為 $s^3 + 4s^2 + 6s +4$
試決定 $G(s)$
Solution:
首先決定誤差轉移函數 $E(s)$ 並將其以 $G(s)$ 與 $R(s)$ 表示：由於 $Y(s) = G(s)E(s)$ 且 $E(s) = R(s) - Y(s)$ 我們可推得
\[
E(s) = R(s) - Y(s) = R(s) - G(s) E(s)
\]故
\[
E(s) = \frac{R(s)}{1 + G(s)}
\]
現在令 $G(s) := \frac{n(s)}{d(s)}$ 則
\[E(s) = \frac{{R(s)}}{{1 + G(s)}} = \frac{{R(s)}}{{1 + \frac{{n\left( s \right)}}{{d\left( s \right)}}}} = \frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}}R(s)
\]故由條件 2 可知
\[
d(s) + n(s) = s^3 + 4s^2 + 6s +4
\]
接著由條件1可知此閉迴路系統對 單位步階訊號 的穩態誤差為零：亦即 $\mathop {\lim }\limits_{s \to 0} sE(s) = 0$，故取 $R(s) = 1/s$ 為單位步階訊號，我們有
\[\mathop {\lim }\limits_{s \to 0} s\frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}}\frac{1}{s} = \mathop {\lim }\limits_{s \to 0} \frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}} = 0\]或者
\[\mathop {\lim }\limits_{s \to 0} \frac{{d\left( s \right)}}{{{s^3} + 4{s^2} + 6s + 4}} = 0
\]且注意到轉移函數必須為真分形式 亦即至少分母階數要大於或等於分子階數，故我們可取
\[
d(s) = s^3 + 4s^2 + 6s
\]作為其中一種選擇。至此我們已經決定 $d(s)$ 又因為 $d(s) + n(s) = s^3 + 4s^2 + 6s + 4$，以上例而言，$n(s) = 4$。故我們得到
\[
G(s) = \frac{n(s)}{d(s)} = \frac{4}{s^3 + 4s^2 + 6s}
\]
現在假設
1. 閉迴路系統對 單位步階訊號 的穩態誤差為零：
2. 閉迴路轉移函數 $Y(s)/R(s)$ 的特性方程式為 $s^3 + 4s^2 + 6s +4$
試決定 $G(s)$
Solution:
首先決定誤差轉移函數 $E(s)$ 並將其以 $G(s)$ 與 $R(s)$ 表示：由於 $Y(s) = G(s)E(s)$ 且 $E(s) = R(s) - Y(s)$ 我們可推得
\[
E(s) = R(s) - Y(s) = R(s) - G(s) E(s)
\]故
\[
E(s) = \frac{R(s)}{1 + G(s)}
\]
現在令 $G(s) := \frac{n(s)}{d(s)}$ 則
\[E(s) = \frac{{R(s)}}{{1 + G(s)}} = \frac{{R(s)}}{{1 + \frac{{n\left( s \right)}}{{d\left( s \right)}}}} = \frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}}R(s)
\]故由條件 2 可知
\[
d(s) + n(s) = s^3 + 4s^2 + 6s +4
\]
接著由條件1可知此閉迴路系統對 單位步階訊號 的穩態誤差為零：亦即 $\mathop {\lim }\limits_{s \to 0} sE(s) = 0$，故取 $R(s) = 1/s$ 為單位步階訊號，我們有
\[\mathop {\lim }\limits_{s \to 0} s\frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}}\frac{1}{s} = \mathop {\lim }\limits_{s \to 0} \frac{{d\left( s \right)}}{{d\left( s \right) + n\left( s \right)}} = 0\]或者
\[\mathop {\lim }\limits_{s \to 0} \frac{{d\left( s \right)}}{{{s^3} + 4{s^2} + 6s + 4}} = 0
\]且注意到轉移函數必須為真分形式 亦即至少分母階數要大於或等於分子階數，故我們可取
\[
d(s) = s^3 + 4s^2 + 6s
\]作為其中一種選擇。至此我們已經決定 $d(s)$ 又因為 $d(s) + n(s) = s^3 + 4s^2 + 6s + 4$，以上例而言，$n(s) = 4$。故我們得到
\[
G(s) = \frac{n(s)}{d(s)} = \frac{4}{s^3 + 4s^2 + 6s}
\]
4/17/2015
[系統理論] 線性系統的 Input/Output to State Stability
回憶在線性系統理論中，我們說 系統為可觀測 (observable) 則可以設計觀察器(observer)。但若今天考慮的是 非線性系統，我們是否能有類似的準則來告訴我們何時可以設計觀察器? 或者說是否能提供非線性系統 類比於系統觀察性( observability) 的條件呢?
答案是肯定的，在非線性系統裡面 我們可用 Incremental Input/output-to-state stability (i-IOSS) 來界定系統是否可觀察。 不過如果是只針對線性系統，則我們僅需要 IOSS 即可等價 observability；以下為定義：
=======================
Definition: Input/output-to-state stability (IOSS)
假設非線性系統 $x^+ = f(x,w); \;\; y=h(x)$ 為 input/output-to-state stable (IOSS) 若下列條件成立：
對任意 $x_0 \in \mathbb{R}^n$ 與 $k \ge 0$，存在 $\beta(\cdot) \in \mathcal{KL}$ 與 $\gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$ 使得
\[
|x(k;x_0, {\bf w})| \le \beta(|x_0|,k) + \gamma_1(||{\bf w}||) + \gamma_2 (||{\bf y}||)
\]其中 $x(k;x_0,{\bf w})$ 為前述非線性系統 在時間 $k$ 與 初始條件 $x_0$ 輸入為 ${\bf w}$的解；另外 $||{\bf w}|| := \max_{j\ge 0} |w(j)|, ||{\bf y}|| := \max_{j\ge 0} |y(j)|$
=======================
那麼有了上述結果，我們應立即想到此定義對於原本線性系統是否有用? 現在我們看個例子：
Example
考慮離散時間 線性系統
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
y = Cx
\end{array}\]試證若系統為 observable，則 系統為 Input-Output State Stable (IOSS)
Proof:
回憶 IOSS 定義：
\[\left| {x\left( {k;{x_0}} \right)} \right| \le \beta \left( {|{x_0}|,k} \right) + {\gamma _1}\left( {||w||} \right) + {\gamma _2}\left( {||y||} \right)
\]其中 $\beta(\cdot) \in \mathcal{KL}$ 且 $\gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$
由於系統為 observable，我們可以設計觀察器，故我們可改寫系統方程如下
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
\Rightarrow {x^ + } = Ax + Gw + Ly - Ly\\
\Rightarrow {x^ + } = \left( {A - LC} \right)x + Gw + Ly
\end{array}\]其中 $ L$ 觀測器的增益矩陣， 且由於系統為 observable，故 $(A-LC)$ 為穩定矩陣；亦即 $eig(A-LC) < 1$ 。現在我們求解 $x(k)$ 如下：
\[\begin{array}{l}
\left| {x\left( k \right)} \right| = \left| {{{\left( {A - LC} \right)}^k}{x_0} + \sum\limits_{j = 0}^{k - 1} {{{\left( {A - LC} \right)}^{k - 1 - j}}\left( {Gw\left( j \right) + Ly\left( j \right)} \right)} } \right|\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le \left| {{{\left( {A - LC} \right)}^k}} \right|\left| {{x_0}} \right| + \sum\limits_{j = 0}^{k - 1} {\left| {{{\left( {A - LC} \right)}^{k - 1 - j}}Gw\left( j \right)} \right|} \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}&{}&{}&{}
\end{array} + \sum\limits_{j = 0}^{k - 1} {\left| {{{\left( {A - LC} \right)}^{k - 1 - j}}Ly\left( j \right)} \right|}
\end{array}\]現在回憶下列結果 (Horn and Johnson, 1985, p.299)：對 $c>0$，
\[\left| {{{\left( {A - LC} \right)}^k}} \right| \le c{\lambda ^k};\begin{array}{*{20}{c}}
{}&{}
\end{array}\mathop {\max }\limits_i \left| {ei{g_i}\left( A - LC \right)} \right| < \lambda < 1
\]故可知
\[\begin{array}{*{20}{l}}
{\left| {x\left( k \right)} \right| \le c{\lambda ^k}\left| {{x_0}} \right| + c{\lambda ^{k - 1}}\sum\limits_{j = 0}^{k - 1} {{\lambda ^{ - j}}\left| {Gw\left( j \right)} \right|} + c{\lambda ^{k - 1}}\sum\limits_{j = 0}^{k - 1} {{\lambda ^{ - j}}\left| {Lw\left( j \right)} \right|} }\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le c{\lambda ^k}\left| {{x_0}} \right| + c\left| G \right|\left\| w \right\|\frac{{1 - {\lambda ^k}}}{{1 - \lambda }} + c\left| L \right|\left\| y \right\|\frac{{1 - {\lambda ^k}}}{{1 - \lambda }}}\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le c{\lambda ^k}\left| {{x_0}} \right| + c\left| G \right|\left\| w \right\|\frac{1}{{1 - \lambda }} + c\left| L \right|\left\| y \right\|\frac{1}{{1 - \lambda }}}
\end{array}\]故我們可取
\[\left\{ \begin{array}{l}
\beta \left( {|{x_0}|,k} \right): = c{\lambda ^k}\left| {{x_0}} \right|\\
{\gamma _1}\left( {||w||} \right): = c\left| G \right|\left\| w \right\|\frac{{ {1 }}}{{1 - \lambda }}\\
{\gamma _2}\left( {||y||} \right){\rm{: = }}c\left| L \right|\left\| y \right\|\frac{{{1 }}}{{1 - \lambda }}
\end{array} \right.\]最後我們僅需檢驗 $\beta(\cdot) \in \mathcal{KL}$, 且 $ \gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$。 $\square$
事實上上述例子可以推廣成如下重要結果：
===================
Theorem: IOSS in linear system is equivalent to detectability
考慮離散時間 線性系統
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
y = Cx
\end{array}\]此系統為 detectable 若且為若 系統為 IOSS
===================
Proof: omitted
答案是肯定的，在非線性系統裡面 我們可用 Incremental Input/output-to-state stability (i-IOSS) 來界定系統是否可觀察。 不過如果是只針對線性系統，則我們僅需要 IOSS 即可等價 observability；以下為定義：
=======================
Definition: Input/output-to-state stability (IOSS)
假設非線性系統 $x^+ = f(x,w); \;\; y=h(x)$ 為 input/output-to-state stable (IOSS) 若下列條件成立：
對任意 $x_0 \in \mathbb{R}^n$ 與 $k \ge 0$，存在 $\beta(\cdot) \in \mathcal{KL}$ 與 $\gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$ 使得
\[
|x(k;x_0, {\bf w})| \le \beta(|x_0|,k) + \gamma_1(||{\bf w}||) + \gamma_2 (||{\bf y}||)
\]其中 $x(k;x_0,{\bf w})$ 為前述非線性系統 在時間 $k$ 與 初始條件 $x_0$ 輸入為 ${\bf w}$的解；另外 $||{\bf w}|| := \max_{j\ge 0} |w(j)|, ||{\bf y}|| := \max_{j\ge 0} |y(j)|$
=======================
那麼有了上述結果，我們應立即想到此定義對於原本線性系統是否有用? 現在我們看個例子：
Example
考慮離散時間 線性系統
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
y = Cx
\end{array}\]試證若系統為 observable，則 系統為 Input-Output State Stable (IOSS)
Proof:
回憶 IOSS 定義：
\[\left| {x\left( {k;{x_0}} \right)} \right| \le \beta \left( {|{x_0}|,k} \right) + {\gamma _1}\left( {||w||} \right) + {\gamma _2}\left( {||y||} \right)
\]其中 $\beta(\cdot) \in \mathcal{KL}$ 且 $\gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$
由於系統為 observable，我們可以設計觀察器，故我們可改寫系統方程如下
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
\Rightarrow {x^ + } = Ax + Gw + Ly - Ly\\
\Rightarrow {x^ + } = \left( {A - LC} \right)x + Gw + Ly
\end{array}\]其中 $ L$ 觀測器的增益矩陣， 且由於系統為 observable，故 $(A-LC)$ 為穩定矩陣；亦即 $eig(A-LC) < 1$ 。現在我們求解 $x(k)$ 如下：
\[\begin{array}{l}
\left| {x\left( k \right)} \right| = \left| {{{\left( {A - LC} \right)}^k}{x_0} + \sum\limits_{j = 0}^{k - 1} {{{\left( {A - LC} \right)}^{k - 1 - j}}\left( {Gw\left( j \right) + Ly\left( j \right)} \right)} } \right|\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le \left| {{{\left( {A - LC} \right)}^k}} \right|\left| {{x_0}} \right| + \sum\limits_{j = 0}^{k - 1} {\left| {{{\left( {A - LC} \right)}^{k - 1 - j}}Gw\left( j \right)} \right|} \\
\begin{array}{*{20}{c}}
{}&{}&{}&{}&{}&{}&{}&{}
\end{array} + \sum\limits_{j = 0}^{k - 1} {\left| {{{\left( {A - LC} \right)}^{k - 1 - j}}Ly\left( j \right)} \right|}
\end{array}\]現在回憶下列結果 (Horn and Johnson, 1985, p.299)：對 $c>0$，
\[\left| {{{\left( {A - LC} \right)}^k}} \right| \le c{\lambda ^k};\begin{array}{*{20}{c}}
{}&{}
\end{array}\mathop {\max }\limits_i \left| {ei{g_i}\left( A - LC \right)} \right| < \lambda < 1
\]故可知
\[\begin{array}{*{20}{l}}
{\left| {x\left( k \right)} \right| \le c{\lambda ^k}\left| {{x_0}} \right| + c{\lambda ^{k - 1}}\sum\limits_{j = 0}^{k - 1} {{\lambda ^{ - j}}\left| {Gw\left( j \right)} \right|} + c{\lambda ^{k - 1}}\sum\limits_{j = 0}^{k - 1} {{\lambda ^{ - j}}\left| {Lw\left( j \right)} \right|} }\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le c{\lambda ^k}\left| {{x_0}} \right| + c\left| G \right|\left\| w \right\|\frac{{1 - {\lambda ^k}}}{{1 - \lambda }} + c\left| L \right|\left\| y \right\|\frac{{1 - {\lambda ^k}}}{{1 - \lambda }}}\\
{\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} \le c{\lambda ^k}\left| {{x_0}} \right| + c\left| G \right|\left\| w \right\|\frac{1}{{1 - \lambda }} + c\left| L \right|\left\| y \right\|\frac{1}{{1 - \lambda }}}
\end{array}\]故我們可取
\[\left\{ \begin{array}{l}
\beta \left( {|{x_0}|,k} \right): = c{\lambda ^k}\left| {{x_0}} \right|\\
{\gamma _1}\left( {||w||} \right): = c\left| G \right|\left\| w \right\|\frac{{ {1 }}}{{1 - \lambda }}\\
{\gamma _2}\left( {||y||} \right){\rm{: = }}c\left| L \right|\left\| y \right\|\frac{{{1 }}}{{1 - \lambda }}
\end{array} \right.\]最後我們僅需檢驗 $\beta(\cdot) \in \mathcal{KL}$, 且 $ \gamma_1(\cdot), \gamma_2(\cdot) \in \mathcal{K}$。 $\square$
事實上上述例子可以推廣成如下重要結果：
===================
Theorem: IOSS in linear system is equivalent to detectability
考慮離散時間 線性系統
\[\begin{array}{l}
{x^ + } = Ax + Gw\\
y = Cx
\end{array}\]此系統為 detectable 若且為若 系統為 IOSS
===================
Proof: omitted
3/14/2015
[系統理論] 離散時間系統的穩定度理論 (1) - Lyapunov Stability Theory
延續前篇 [系統理論] 離散時間系統的穩定度理論 (0) - 先備概念，我們現在可以開始介紹 Lyapunov Stability Theory。
Definition: Lyapunov Function
一個函數 $V: \mathbb{R}^n \to \mathbb{R}_{\ge 0}$ 被稱作為 Lyapunov function for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$ 若下列條件成立：
存在函數 $\alpha_1(\cdot), \alpha_2(\cdot), \alpha_3(\cdot) \in \mathcal{K}_\infty$ 使得對任意 $x \in \mathbb{R}^n$，
===================
Definition: Lyapunov Function
一個函數 $V: \mathbb{R}^n \to \mathbb{R}_{\ge 0}$ 被稱作為 Lyapunov function for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$ 若下列條件成立：
存在函數 $\alpha_1(\cdot), \alpha_2(\cdot), \alpha_3(\cdot) \in \mathcal{K}_\infty$ 使得對任意 $x \in \mathbb{R}^n$，
- $V(x) \ge \alpha_1(|x|_\mathcal{A})$
- $V(x) \le \alpha_2(|x|_\mathcal{A})$
- $V(f(x)) - V(x) \le -\alpha_3(|x|_\mathcal{A})$
給定 $\mathcal{A}$ 為 closed positive invariant for $x^+ = f(x)$ 且 $\mathcal{A} \subset X$ 我們說函數 $V(\cdot)$ 為 Lyapunov function in $X$ for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$ 若下列條件成立：
對任意 $x \in X$，$V(\cdot)$ 滿足上述三條不等式。
上述 Lyapunov function 與 globally asymptotically stable 息息相關，事實上此Lyapunov function 的存在性為 globally asymptotically stable 的充分條件。我們將此記做下方結果
===================
Theorem 1: Existence of Lyapunov Function Implies Globally Asymptotically Stability
假設 $V(\cdot)$ 為 Lyapuonv function for $x^+ = f(x)$ 與 $\mathcal{A}$，則 $\mathcal{A}$ 為 globally asymptotically stable。
===================
Proof:
我們要證 $\mathcal{A}$ 為 globally asymptotically stable。故須證明
1. $\mathcal{A}$ 為 locally stable。
2. $\mathcal{A}$ 為 globally attractive。
先證 $\mathcal{A}$ 為 locally stable：給定 $\varepsilon >0$ 我們要找出 $\delta >0$ 使得對任意 $i \in \mathbb{Z}_{\ge 0}$ $|x|_\mathcal{A}< \delta \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $
由於 $V$ 為 Lyapunov function 故由其定義可繪製下圖幫助我們選擇 $\delta$
由於 $V$ 為 Lyapunov function 故由其定義可繪製下圖幫助我們選擇 $\delta$
令 $\delta : = \alpha _2^{ - 1}\left( {{\alpha _1}\left( \varepsilon \right)} \right)$；現在給定 $i \in \mathbb{Z}_{\ge 0}$， 我們要證明 $ |x|_\mathcal{A}< \delta \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $ ；故假設 $|x|_\mathcal{A}< \delta$ 則結合前述 $\delta$ 定義 我們有
\[{\left| x \right|_{{\cal A}}} < \delta = \alpha _2^{ - 1}\left( {{\alpha _1}\left( \varepsilon \right)} \right)\]故
\[\Rightarrow \alpha _2^{}\left( {{{\left| x \right|}_{{\cal A}}}} \right) < {\alpha _1}\left( \varepsilon \right)\]由Lyapunov function 定義第2條不等式可知
\[V(x) \le {\alpha _2}(|x{|_{{\cal A}}}) \Rightarrow V(x) \le {\alpha _1}\left( \varepsilon \right) \ \ \ \ (*)
\]現在觀察 Lyapunov function 定義第3條不等式，並且令 $\phi(i,x)$ 為 $x^+ = f(x)$ 之解，且注意到 $\alpha_i(\cdot)$ 其中 $i=1,2,3$ 皆為 $\mathcal{K}_\infty$ 函數，故我們可寫下
\[\begin{array}{*{20}{l}}
{V(f(x)) - V(x) \le - {\alpha _3}\left( {|x|} \right)}\\
{ \Rightarrow V(f(x)) \le V(x)}\\
{ \Rightarrow V(f(x)) \le V(x) \le {\alpha _1}\left( \varepsilon \right)\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}by\begin{array}{*{20}{c}}
{}
\end{array}\left( * \right).}\\
{ \Rightarrow V(\phi \left( {i;x} \right)) \le {\alpha _1}\left( \varepsilon \right)}\\
{ \Rightarrow {\alpha _1}(|\phi \left( {i;x} \right){|_A}) \le V(\phi \left( {i;x} \right)) \le {\alpha _1}\left( \varepsilon \right)\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}by\begin{array}{*{20}{c}}
{}
\end{array}1.}\\
{ \Rightarrow |\phi \left( {i;x} \right){|_A} \le \varepsilon }
\end{array}\]至此我們得到 $|x|_\mathcal{A}< \delta \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $ 即為所求。
接著我們證明 global attractivity：故給定任意 $x \in X$， 要證明 $|\phi(i;x)|_\mathcal{A} \to 0 \text{ as $i \to \infty$}$。基本想法為比較不同時間點 $\phi$。 現在給定 $\phi(i;x)$ 為 $x^+ = f(x)$ 之解，由 Lyapunov function 定義第3條不等式可知
\[\begin{array}{l}
V(f(x)) - V(x) \le - {\alpha _3}\left( {\left| x \right|} \right)\\
\Rightarrow V(\phi \left( {i + 1;x} \right)) - V(\phi \left( {i;x} \right)) \le - {\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)
\end{array}\]令 $V_{i+1}:= V(\phi \left( {i + 1;x} \right))$ 且 $V_i := V(\phi \left( {i;x} \right))$ 則由上式可推知 對任意 $x$ 而言， 數列 $\{V\}_i$ 為非遞增數列 且 有下界為 $0$，故可知 $V_i$ 收斂亦即
\[{V_{i + 1}} - {V_i} \le - {\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right) \to 0\begin{array}{*{20}{c}}
{}
\end{array}as\begin{array}{*{20}{c}}
{}
\end{array}i \to \infty \]故${\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right) \to 0$ 又由於 $\alpha_3 \in \mathcal{K}_\infty$，故
\[\begin{array}{l}
\left| {\phi \left( {i;x} \right)} \right| = \alpha _3^{ - 1}\underbrace {\left( {{\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)} \right)}_{ \to 0}\\
\Rightarrow \left| {\phi \left( {i;x} \right)} \right| = \alpha _3^{ - 1}\left( {{\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)} \right) \to 0\begin{array}{*{20}{c}}
{}&{}&{}
\end{array}since\begin{array}{*{20}{c}}
{}
\end{array}\alpha _3^{ - 1} \in {\mathcal{K}_\infty }
\end{array}\]至此證畢。 $\square$
上述定理告訴我們 asymptotically stable 的充分條件，但對於 必要條件 並無著墨。所幸透過適度的增強假設，我們仍可得到 asymptotically stable 必要條件，在此紀錄如下：
=============
=============
故我們將 Theorem 1 與 Theorem 2 整合可得如下充分必要條件：
=============
Theorem 3
若 $\mathcal{A}$ 為 compact 且 $f(\cdot)$ 為連續函數。則集合 $\mathcal{A}$ 為 globally asymptotic stable for the system $x^+ = f(x)$ 若且唯若 存在 平滑 (smooth) Lyapunov function for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$。
=============
\[{\left| x \right|_{{\cal A}}} < \delta = \alpha _2^{ - 1}\left( {{\alpha _1}\left( \varepsilon \right)} \right)\]故
\[\Rightarrow \alpha _2^{}\left( {{{\left| x \right|}_{{\cal A}}}} \right) < {\alpha _1}\left( \varepsilon \right)\]由Lyapunov function 定義第2條不等式可知
\[V(x) \le {\alpha _2}(|x{|_{{\cal A}}}) \Rightarrow V(x) \le {\alpha _1}\left( \varepsilon \right) \ \ \ \ (*)
\]現在觀察 Lyapunov function 定義第3條不等式，並且令 $\phi(i,x)$ 為 $x^+ = f(x)$ 之解，且注意到 $\alpha_i(\cdot)$ 其中 $i=1,2,3$ 皆為 $\mathcal{K}_\infty$ 函數，故我們可寫下
\[\begin{array}{*{20}{l}}
{V(f(x)) - V(x) \le - {\alpha _3}\left( {|x|} \right)}\\
{ \Rightarrow V(f(x)) \le V(x)}\\
{ \Rightarrow V(f(x)) \le V(x) \le {\alpha _1}\left( \varepsilon \right)\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}by\begin{array}{*{20}{c}}
{}
\end{array}\left( * \right).}\\
{ \Rightarrow V(\phi \left( {i;x} \right)) \le {\alpha _1}\left( \varepsilon \right)}\\
{ \Rightarrow {\alpha _1}(|\phi \left( {i;x} \right){|_A}) \le V(\phi \left( {i;x} \right)) \le {\alpha _1}\left( \varepsilon \right)\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array}by\begin{array}{*{20}{c}}
{}
\end{array}1.}\\
{ \Rightarrow |\phi \left( {i;x} \right){|_A} \le \varepsilon }
\end{array}\]至此我們得到 $|x|_\mathcal{A}< \delta \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $ 即為所求。
接著我們證明 global attractivity：故給定任意 $x \in X$， 要證明 $|\phi(i;x)|_\mathcal{A} \to 0 \text{ as $i \to \infty$}$。基本想法為比較不同時間點 $\phi$。 現在給定 $\phi(i;x)$ 為 $x^+ = f(x)$ 之解，由 Lyapunov function 定義第3條不等式可知
\[\begin{array}{l}
V(f(x)) - V(x) \le - {\alpha _3}\left( {\left| x \right|} \right)\\
\Rightarrow V(\phi \left( {i + 1;x} \right)) - V(\phi \left( {i;x} \right)) \le - {\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)
\end{array}\]令 $V_{i+1}:= V(\phi \left( {i + 1;x} \right))$ 且 $V_i := V(\phi \left( {i;x} \right))$ 則由上式可推知 對任意 $x$ 而言， 數列 $\{V\}_i$ 為非遞增數列 且 有下界為 $0$，故可知 $V_i$ 收斂亦即
\[{V_{i + 1}} - {V_i} \le - {\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right) \to 0\begin{array}{*{20}{c}}
{}
\end{array}as\begin{array}{*{20}{c}}
{}
\end{array}i \to \infty \]故${\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right) \to 0$ 又由於 $\alpha_3 \in \mathcal{K}_\infty$，故
\[\begin{array}{l}
\left| {\phi \left( {i;x} \right)} \right| = \alpha _3^{ - 1}\underbrace {\left( {{\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)} \right)}_{ \to 0}\\
\Rightarrow \left| {\phi \left( {i;x} \right)} \right| = \alpha _3^{ - 1}\left( {{\alpha _3}\left( {\left| {\phi \left( {i;x} \right)} \right|} \right)} \right) \to 0\begin{array}{*{20}{c}}
{}&{}&{}
\end{array}since\begin{array}{*{20}{c}}
{}
\end{array}\alpha _3^{ - 1} \in {\mathcal{K}_\infty }
\end{array}\]至此證畢。 $\square$
上述定理告訴我們 asymptotically stable 的充分條件，但對於 必要條件 並無著墨。所幸透過適度的增強假設，我們仍可得到 asymptotically stable 必要條件，在此紀錄如下：
=============
Theorem 2: Converse Theorem for Asymptotic Stability
令 $\mathcal{A}$ 為 compact 且 $f(\cdot)$ 為連續函數。假設 $\mathcal{A}$ 為 globally asymptotic stable for the system $x^+ = f(x)$ 則 存在 平滑 (smooth) Lyapunov function for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$。=============
Proof: omitted
=============
Theorem 3
若 $\mathcal{A}$ 為 compact 且 $f(\cdot)$ 為連續函數。則集合 $\mathcal{A}$ 為 globally asymptotic stable for the system $x^+ = f(x)$ 若且唯若 存在 平滑 (smooth) Lyapunov function for system $x^+ = f(x)$ 與 集合 $\mathcal{A}$。
=============
[系統理論] 離散時間系統的穩定度理論 (0) - 先備概念
穩定度理論可追朔至 Aleksandr Lyapunov在 1892 出版的 The General Problem of Stability of Motion 提出，主要是透過建構 Lyapunov 函數 來判別動態系統是否穩定。以下討論我們將以 離散時間 非線性動態系統為主。
現在考慮以下 離散時間 非線性動態系統
\[
x^+ = f(x,u)
\]其中 $x \in \mathbb{R}^n$ 為當前系統狀態 且 $u \in \mathbb{R}^m$ 為當前的控制力；$x^+$ 為下個時刻的系統狀態。且假設 $f: \mathbb{R}^n \times \mathbb{R}^m \to \mathbb{R}^n$ 為連續函數。
定義 $\phi(k; x, {\bf u}) $為在時刻 $k$，對於動態系統 $x^+ = f(x,u)$ 的解 (初始值為 $x(0)=x$ ； 控制力序列 $\bf u$ $:=\{u(0), u(1), ...\}$)
若 控制律 $u := \kappa (x)$ 決定，則系統閉迴路可表為
\[
x^+ = f(x,\kappa(x)):=f_c(x)
\]注意到 $\kappa(\cdot)$ 不一定為連續函數，此時對應的 $f(x, \kappa(\cdot))$ 亦不一定為連續。對此不連續的情況我們額外假設 $f_c(\cdot)$ 為 局部有界(locally bounded)。
目標：我們希望 控制系統 要"穩定"。
在此所謂的穩定 意指 控制系統對於 初始狀態 的小擾動 不會 導致 閉迴路系統響應 大幅度擾動 且 系統狀態能夠收斂到指定的狀態 或者 收斂到指定的 狀態集合 (此情況多半發生在有外部干擾的時候)。
以下我們會針對定義 系統的 穩定度 與 漸進穩定度；在介紹之前我們需要先定義一些名詞：首先是 如何指出系統狀態的收斂
============
Definition: Equilibrium Point or Steady-State
狀態 $x^*$ 被稱作 $x^+ = f(x)$ 的 平衡點(equilibrium point) 若 $$x(0) = x^* \Rightarrow x(k) = \phi(k;x^*) = x^*, \;\; \forall k \ge 0$$
============
Comment:
1. 上述定義表示 $x^*$ 為 平衡點 若其滿足 $x^* = f(x^*)$
2. equilibrium point 為 被隔離的(isolated) 若在 $x^*$ 附近沒有其他的平衡點。
3. 非線性系統可能有多個 被隔離的平衡點
Example: Equilibrium Point of Linear System
考慮離散時間線性系統
\[
x^+ = Ax +b
\]具有平衡點 $x^*$ 則
\[\begin{array}{l}
{x^ + } = Ax + b\
\Rightarrow {x^*} = A{x^*} + b\\
\Rightarrow \left( {I - A} \right){x^*} = b
\end{array}\]若 $I-A$ 反矩陣存在，則我們說 此線性系統有 unique (isolated) 平衡點
\[{x^*} = {\left( {I - A} \right)^{ - 1}}b\]若 $I-A$ 反矩陣不存在，則我們說此線性系統有 連續統 (continuum) $\{x: (I-A) x=b\}$ 的平衡點。
若我們考慮 震盪系統 的穩定度，則此時不再是討論 是否收斂到某個狀態 (平衡點)；而是討論收斂到某個集合。以下我們給出此類集合所需的定義：
一個集合 $\mathcal{A}$ 稱作 positive invariant for system $x^+ = f(x)$ 若下列條件成立：
\[
x \in A \Rightarrow f(x) \in \mathcal{A}
\]=============
一個函數 $g: \mathbb{R}_{\ge 0} \to \mathbb{R}_{\ge 0}$ 為 $\mathcal{K}$ 類函數若下列條件滿足：
我們說 一個函數 $g: \mathbb{R}_{\ge 0} \to \mathbb{R}_{\ge 0}$ 為 $\mathcal{K}_\infty$ 類函數若下列條件滿足：
FACT 2:
若 $\alpha_1(\cdot)$ 與 $\alpha_2(\cdot)$ 為 $\mathcal{K}$ 類函數 且 $\beta(\cdot)$ 為 $\mathcal{KL}$ 函數，則 $\sigma (r,s): = {\alpha _1}(\beta \left( {{\alpha _2}\left( r \right)} \right),s)$ 為 $\mathcal{KL}$函數。
現在考慮以下 離散時間 非線性動態系統
\[
x^+ = f(x,u)
\]其中 $x \in \mathbb{R}^n$ 為當前系統狀態 且 $u \in \mathbb{R}^m$ 為當前的控制力；$x^+$ 為下個時刻的系統狀態。且假設 $f: \mathbb{R}^n \times \mathbb{R}^m \to \mathbb{R}^n$ 為連續函數。
定義 $\phi(k; x, {\bf u}) $為在時刻 $k$，對於動態系統 $x^+ = f(x,u)$ 的解 (初始值為 $x(0)=x$ ； 控制力序列 $\bf u$ $:=\{u(0), u(1), ...\}$)
若 控制律 $u := \kappa (x)$ 決定，則系統閉迴路可表為
\[
x^+ = f(x,\kappa(x)):=f_c(x)
\]注意到 $\kappa(\cdot)$ 不一定為連續函數，此時對應的 $f(x, \kappa(\cdot))$ 亦不一定為連續。對此不連續的情況我們額外假設 $f_c(\cdot)$ 為 局部有界(locally bounded)。
目標：我們希望 控制系統 要"穩定"。
在此所謂的穩定 意指 控制系統對於 初始狀態 的小擾動 不會 導致 閉迴路系統響應 大幅度擾動 且 系統狀態能夠收斂到指定的狀態 或者 收斂到指定的 狀態集合 (此情況多半發生在有外部干擾的時候)。
以下我們會針對定義 系統的 穩定度 與 漸進穩定度；在介紹之前我們需要先定義一些名詞：首先是 如何指出系統狀態的收斂
============
Definition: Equilibrium Point or Steady-State
狀態 $x^*$ 被稱作 $x^+ = f(x)$ 的 平衡點(equilibrium point) 若 $$x(0) = x^* \Rightarrow x(k) = \phi(k;x^*) = x^*, \;\; \forall k \ge 0$$
============
1. 上述定義表示 $x^*$ 為 平衡點 若其滿足 $x^* = f(x^*)$
2. equilibrium point 為 被隔離的(isolated) 若在 $x^*$ 附近沒有其他的平衡點。
3. 非線性系統可能有多個 被隔離的平衡點
Example: Equilibrium Point of Linear System
考慮離散時間線性系統
\[
x^+ = Ax +b
\]具有平衡點 $x^*$ 則
\[\begin{array}{l}
{x^ + } = Ax + b\
\Rightarrow {x^*} = A{x^*} + b\\
\Rightarrow \left( {I - A} \right){x^*} = b
\end{array}\]若 $I-A$ 反矩陣存在，則我們說 此線性系統有 unique (isolated) 平衡點
\[{x^*} = {\left( {I - A} \right)^{ - 1}}b\]若 $I-A$ 反矩陣不存在，則我們說此線性系統有 連續統 (continuum) $\{x: (I-A) x=b\}$ 的平衡點。
若我們考慮 震盪系統 的穩定度，則此時不再是討論 是否收斂到某個狀態 (平衡點)；而是討論收斂到某個集合。以下我們給出此類集合所需的定義：
=============
Definition: Positive Invariant Set一個集合 $\mathcal{A}$ 稱作 positive invariant for system $x^+ = f(x)$ 若下列條件成立：
\[
x \in A \Rightarrow f(x) \in \mathcal{A}
\]=============
Comment:
1. Positive 來自於 $x^+ = f(x)$ 為動態系統隨時間 $k$ "增加" 而持續變動。
2. 考慮 closed set $\mathcal{A}:=\{x^*\}$ 且 $x^*$ 為系統 $x^+ = f(x)$ 的平衡點，則
2. 考慮 closed set $\mathcal{A}:=\{x^*\}$ 且 $x^*$ 為系統 $x^+ = f(x)$ 的平衡點，則
\[
x \in \mathcal{A}\;\; (\text{since} \;x^* \in \mathcal{A}) \Rightarrow f(x) \in \mathcal{A} \;\;(\text{since} \;f(x) = x^*)
\]
x \in \mathcal{A}\;\; (\text{since} \;x^* \in \mathcal{A}) \Rightarrow f(x) \in \mathcal{A} \;\;(\text{since} \;f(x) = x^*)
\]
=============
Definition: K, K infinity, KL function一個函數 $g: \mathbb{R}_{\ge 0} \to \mathbb{R}_{\ge 0}$ 為 $\mathcal{K}$ 類函數若下列條件滿足：
- $g$ 為連續
- $g(0) = 0$
- 嚴格遞增(strictly increasing)；亦即 $\forall x,y$，$y > x \Rightarrow g(y) > g(x)$
我們說 一個函數 $g: \mathbb{R}_{\ge 0} \to \mathbb{R}_{\ge 0}$ 為 $\mathcal{K}_\infty$ 類函數若下列條件滿足：
- $g$ 為 $\mathcal{K}$ 類函數
- 當 $t \to \infty$，$g(t) \to \infty$
我們說一個函數 $h: \mathbb{R}_{\ge 0} \times \mathbb{Z}_{\ge 0} \to \mathbb{R}_{\ge 0}$ 為 $\mathcal{KL}$ 類函數若下列條件滿足：
- 對任意 $t \ge 0$，$h(\cdot, t)$ 為 $\mathcal{K}$ 類函數
- 對任意 $s \ge 0$，$h(s, \cdot)$ 為非遞增(nonincreasing) 且 滿足 $\lim_{t\to \infty} h(s,t) =0$
=============
Example
1. $g(x) := x$ 為 $\mathcal{K}$類函數 (亦為 $\mathcal{K}_\infty $ 函數)
2. $erf(x)$ 為 $\mathcal{K}$類函數
以下我們將前述 $\mathcal{K}$類函數的重要性質：
=============
FACT 1: Inverse K function is a K function
若 $\alpha_1(\cdot), \alpha_2(\cdot)$ 為 $\mathcal{K}$ 類函數 (或者 $\mathcal{K}_\infty$ 函數)，則其反函數 $\alpha_1^{-1}(\cdot), \alpha_1^{-1}(\cdot)$ 亦仍為 $\mathcal{K}$ 類函數 (或者 $\mathcal{K}_\infty$ 函數)
若 $\alpha_1(\cdot), \alpha_2(\cdot)$ 為 $\mathcal{K}$ 類函數 (或者 $\mathcal{K}_\infty$ 函數)，則其反函數 $\alpha_1^{-1}(\cdot), \alpha_1^{-1}(\cdot)$ 亦仍為 $\mathcal{K}$ 類函數 (或者 $\mathcal{K}_\infty$ 函數)
=============
Proof: omitted
=============
若 $\alpha_1(\cdot)$ 與 $\alpha_2(\cdot)$ 為 $\mathcal{K}$ 類函數 且 $\beta(\cdot)$ 為 $\mathcal{KL}$ 函數，則 $\sigma (r,s): = {\alpha _1}(\beta \left( {{\alpha _2}\left( r \right)} \right),s)$ 為 $\mathcal{KL}$函數。
=============
Proof: omitted
有了以上定義我們可以開始引入 穩定度 的嚴格定義。以下我們考慮 $x^+ = f(x)$ 且假設 $f(\cdot)$ 為 局部有界(locally bounded) 且集合 $A$ 為 closed 與 positive invariant 。
==================
Definition: Local Stability (Stability in Lyapunov Sense)
給定 closed positive invariant 集合 $\mathcal{A}$ 。我們稱 此集合 $\mathcal{A}$ 為 locally stable for $x^+ = f(x)$ 若下列條件成立：
對任意 $\varepsilon>0$ 存在 $\delta >0$ 使得對任意 $i \in Z_{\ge 0}$， $|x|_\mathcal{A} < \delta \Rightarrow |\phi(i; x)|_\mathcal{A} < \varepsilon$
其中 $|x|_\mathcal{A} := \inf_{z \in \mathcal{A}} |x - z|$
有了以上定義我們可以開始引入 穩定度 的嚴格定義。以下我們考慮 $x^+ = f(x)$ 且假設 $f(\cdot)$ 為 局部有界(locally bounded) 且集合 $A$ 為 closed 與 positive invariant 。
==================
Definition: Local Stability (Stability in Lyapunov Sense)
給定 closed positive invariant 集合 $\mathcal{A}$ 。我們稱 此集合 $\mathcal{A}$ 為 locally stable for $x^+ = f(x)$ 若下列條件成立：
對任意 $\varepsilon>0$ 存在 $\delta >0$ 使得對任意 $i \in Z_{\ge 0}$， $|x|_\mathcal{A} < \delta \Rightarrow |\phi(i; x)|_\mathcal{A} < \varepsilon$
其中 $|x|_\mathcal{A} := \inf_{z \in \mathcal{A}} |x - z|$
==================
Example:
考慮 $A:= \{0\}$ 則 Local Stability 可由下圖得知
上圖顯示了若給定任意初始位置 $x$ 且此 $x$ 與原點 $A:=\{0\}$ 距離落在 開球 $B_{\delta}$ 之中，且若系統 $x^+ = f(x)$ 隨時間變化演進，其解 $\phi(i,x)$ 到原點距離 $A=\{0\}$ 持續落在另一開球 $B_\varepsilon$之中，故此系統稱為 Local stable。
======================
Definition: Global Attraction
給定 closed positive invariant 集合 $\mathcal{A}$ 。我們說此集合 $A$ 為 globally attractive for system $x^+ = f(x)$ 若下列條件成立：
\[
|\phi(i;x)|_\mathcal{A} \to 0 \text{ as $i \to \infty$} \;\; \forall x\in \mathbb{R}^n
\]======================
======================
Comment:
考慮 $\mathcal{A}:=\{0\}$，有可能 globally attractive 但並非 locally stable。比如說考慮\[
x^+ = Ax + \phi(x)
\]其中 $A$ 有 eigenvalue $\lambda_1 = 0.5$ 與 $\lambda_2 = 2$ 且對應的 eigenvector 為 $w_1, w_2$且 $\phi(\cdot)$為 平滑函數 滿足 $\phi(0) = 0$ 與 ${\left. {\frac{\partial }{{\partial x}}\phi (x)} \right|_{x = 0}} = 0$
故在 $0$ 附近，$x^+ = Ax + \phi(x)$ 行為將會非常接近 $x^+ = Ax$ ；故若 $\phi(x) =0$ 則 特徵向量$w_1$ 因為具有特徵值 $\lambda_1 = 0.5$ (落在 unit circle 之中)故此特徵向量會迫使狀態收斂到 $0$點，但 特徵向量 $w_2$ 具有不穩定的特徵值 $\lambda_2 = 2$ 故此特徵向量會迫使狀態發散。故總和此兩者，可知儘管有 globally attractive 但卻沒有 stable origin ($\mathcal{A}:=\{0\}$無法滿足 local stability 定義)
以下我們將相關的穩定度定義總結如下：
=============================
Definition: Stability without constraint
給定 closed positive invariant 集合 $\mathcal{A}$ 為
以下結果將前述穩定度定義 與 KL 函數做連結：
=============================
FACT: (Globally Asymptotic Stable and KL function)
令集合 $\mathcal{A}$ 為 compact 且 positive invariant； $f(\cdot)$ 為連續函數。則 $\mathcal{A}$ 為 globally asymptotic stable for $x^+ = f(x)$ 若且唯若 存在 $\mathcal{KL}$ 函數 $\beta(\cdot)$ 使得 對任意 $x \in \mathbb{R}^n$，
\[
|\phi(i;x)|_\mathcal{A} \le \beta(|x|_\mathcal{A},i)\;\; \forall i \in \mathbb{Z}_{\ge 0}
\]=============================
另外，實際上若考慮系統狀態有拘束的情形，則 globally asymptotic stability 並不保證能夠達成，此時我們需要再次拓展前述定義來滿足有拘束的情況：
=============================
Definition: Stability with Constraint Set X
假設狀態拘束集合 $X \subset \mathbb{R}^n$ 為 positive invariant for $x^+ = f(x)$ 且集合 $\mathcal{A}$ closed positive invariant for $x^+ = f(x)$ 且 $\mathcal{A} \subset int(X)$ ( $int(X) :=$ interior of $X$) 則我們說集合 $\mathcal{A}$ 為
延伸閱讀
[系統理論] 離散時間系統的穩定度理論 (1) - Lyapunov Stability Theory
ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design", 2009
Example:
考慮 $A:= \{0\}$ 則 Local Stability 可由下圖得知
======================
Definition: Global Attraction
給定 closed positive invariant 集合 $\mathcal{A}$ 。我們說此集合 $A$ 為 globally attractive for system $x^+ = f(x)$ 若下列條件成立：
\[
|\phi(i;x)|_\mathcal{A} \to 0 \text{ as $i \to \infty$} \;\; \forall x\in \mathbb{R}^n
\]======================
======================
Definition: (Global Asymptotic Stability (GAS))
給定 closed positive invariant 集合 $\mathcal{A}$ 為 globally asymptotically stable for system $x^+ = f(x)$ 若下列條件成立：
$\mathcal{A}$ 為 locally stable 且 globally attractive======================
Comment:
考慮 $\mathcal{A}:=\{0\}$，有可能 globally attractive 但並非 locally stable。比如說考慮\[
x^+ = Ax + \phi(x)
\]其中 $A$ 有 eigenvalue $\lambda_1 = 0.5$ 與 $\lambda_2 = 2$ 且對應的 eigenvector 為 $w_1, w_2$且 $\phi(\cdot)$為 平滑函數 滿足 $\phi(0) = 0$ 與 ${\left. {\frac{\partial }{{\partial x}}\phi (x)} \right|_{x = 0}} = 0$
故在 $0$ 附近，$x^+ = Ax + \phi(x)$ 行為將會非常接近 $x^+ = Ax$ ；故若 $\phi(x) =0$ 則 特徵向量$w_1$ 因為具有特徵值 $\lambda_1 = 0.5$ (落在 unit circle 之中)故此特徵向量會迫使狀態收斂到 $0$點，但 特徵向量 $w_2$ 具有不穩定的特徵值 $\lambda_2 = 2$ 故此特徵向量會迫使狀態發散。故總和此兩者，可知儘管有 globally attractive 但卻沒有 stable origin ($\mathcal{A}:=\{0\}$無法滿足 local stability 定義)
以下我們將相關的穩定度定義總結如下：
=============================
Definition: Stability without constraint
給定 closed positive invariant 集合 $\mathcal{A}$ 為
- locally stable 若 對任意 $\varepsilon >0$ 存在 $\delta >0$ 使得對任意 $i \in \mathbb{Z}_{\ge 0}$ $|x|_\mathcal{A} < \delta \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $
- unstable 若 其 不為 locally stable
- locally attractive 若 存在 $\eta >0$ 使得 $|x|_\mathcal{A} < \eta \Rightarrow |\phi(i;x)|_\mathcal{A} \to 0\;\; \text{as $i \to \infty$} $
- globally attractive 若對任意 $x \in \mathbb{R}^n$， $|\phi(i;x)|_\mathcal{A} \to 0\;\; \text{as $i \to \infty$} $
- locally asymptotically stable 若其為 locally stable 與 locally attractive
- globally asymptotically stable 若其為 locally stable 與 globally attractive
- locally exponentially stable 若 存在 $\eta >0, c>0$ 與 $\gamma \in (0,1)$ 使得 對任意 $i \in \mathbb{Z}_{\ge 0}$ 而言， $|x|_\mathcal{A} < \eta \Rightarrow |\phi(i;x)|_{\mathcal{A}} \le c |x|_\mathcal{A} \gamma^i$
- globally exponentially stable 若 存在 $c>0$ 與 $\gamma \in (0,1)$ 使得 對任意 $i \in \mathbb{Z}_{\ge 0}$ 而言， $|\phi(i;x)|_{\mathcal{A}} \le c |x|_\mathcal{A} \gamma^i$
=============================
FACT: (Globally Asymptotic Stable and KL function)
令集合 $\mathcal{A}$ 為 compact 且 positive invariant； $f(\cdot)$ 為連續函數。則 $\mathcal{A}$ 為 globally asymptotic stable for $x^+ = f(x)$ 若且唯若 存在 $\mathcal{KL}$ 函數 $\beta(\cdot)$ 使得 對任意 $x \in \mathbb{R}^n$，
\[
|\phi(i;x)|_\mathcal{A} \le \beta(|x|_\mathcal{A},i)\;\; \forall i \in \mathbb{Z}_{\ge 0}
\]=============================
另外，實際上若考慮系統狀態有拘束的情形，則 globally asymptotic stability 並不保證能夠達成，此時我們需要再次拓展前述定義來滿足有拘束的情況：
=============================
Definition: Stability with Constraint Set X
假設狀態拘束集合 $X \subset \mathbb{R}^n$ 為 positive invariant for $x^+ = f(x)$ 且集合 $\mathcal{A}$ closed positive invariant for $x^+ = f(x)$ 且 $\mathcal{A} \subset int(X)$ ( $int(X) :=$ interior of $X$) 則我們說集合 $\mathcal{A}$ 為
- locally stable in $X$ 若 對任意 $\varepsilon >0$ 存在 $\delta >0$ 使得對任意 $i \in \mathbb{Z}_{\ge 0}$ $x \in X \cap (\mathcal{A} \oplus B_\delta) \Rightarrow |\phi(i;x)|_\mathcal{A} <\varepsilon $
- locally attractive in $X$ 若 存在 $\eta >0$ 使得 $x \in X \cap (\mathcal{A} \oplus B_\delta) \Rightarrow |\phi(i;x)|_\mathcal{A} \to 0\;\; \text{as $i \to \infty$} $
- attractive in $X$ 若 $ |\phi(i;x)|_\mathcal{A} \to 0\;\; \text{as $i \to \infty$} \; \forall x \in X$
- locally asymptotically stable in $X$ 若其為 locally stable in $X$ 與 locally attractive in $X$
- asymptotically stable with region of attraction $X$ 若其為 locally stable in $X$ 與 attractive in $X$。
- locally exponentially stable with region of attraction $X$ 若 存在 $\eta >0, c>0$ 與 $\gamma \in (0,1)$ 使得 對任意 $i \in \mathbb{Z}_{\ge 0}$ 而言， $x \in X \cap (\mathcal{A} \oplus B_\eta) \Rightarrow |\phi(i;x)|_{\mathcal{A}} \le c |x|_\mathcal{A} \gamma^i$
- globally exponentially stable with region of attraction $X$ 若 存在 $c>0$ 與 $\gamma \in (0,1)$ 使得 對任意 $i \in \mathbb{Z}_{\ge 0}$ 而言， $|\phi(i;x)|_{\mathcal{A}} \le c |x|_\mathcal{A} \gamma^i$
[系統理論] 離散時間系統的穩定度理論 (1) - Lyapunov Stability Theory
ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design", 2009
2/21/2015
[最佳估計] 狀態估測器的收斂性
考慮離散時間 LTI 系統
\[\begin{array}{l}
x(k + 1) = Ax(k) + w(k)\\
y(k) = Cx\left( k \right) + v(k)
\end{array}\]其中 $x \in \mathbb{R}^n, u \in \mathbb{R}^m, y \in \mathbb{R}^p$ 且 $x(0) = x_0$；
基本狀態估計問題：給定初始估計誤差 ( i.e., $\bar x(0) \neq x_0$) 且 不考慮 雜訊 無外部干擾的情況 ($w(k) = v(k)=0$)，我們想問 $\hat x(k) \to x(k)$ as $k \to \infty$ ?
=================
Theorem: (Convergence of Estimator Cost)
給定無 noise 量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$ 則最佳估測器的 cost $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 。
=================
Proof:
由於
\[V_T^{} = \frac{1}{2}\left( \begin{array}{l}
\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{k = 0}^{T - 1} {|\hat x\left( {k + 1} \right) - A\hat x\left( k \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{k = 0}^T {\left| {y\left( k \right) - C\hat x\left( k \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right)\]故我們分三個步驟證明 $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 ：
首先證明 sequence $\{V_T^*\}$ 有上界(bounded above)；接著我們證明 $\{V_T^*\}$ 為 nondecreasing。則由步驟一與步驟二可推論 $\{V_T^*\}$ 必定收斂。
現在我們開始證明 $V_T^*({\bf y}(T)) $ 有界：
由於我們的目標為 $\min_{\hat{ {\bf x}}(k)} V_k$ 但我們並知道 怎樣的 $\hat x$ 可以幫助我們達成此目標，故我們先暫取 $\hat x(0) := x_0$ (並不一定為最佳解!) 則我們有
\[\begin{array}{l}
\left\{ \begin{array}{l}
\hat x(1) = A\hat x(0)\\
\hat x(2) = A\hat x(1) = {A^2}\hat x(0)\\
\vdots
\end{array} \right.\\
\Rightarrow y\left( k \right) = C{A^k}{x_0}
\end{array}\]故若將此結果帶入我們的 cost 可得
\[\bar{V}_T^{} = \frac{1}{2}\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2 < \infty
\]但注意到此並非最佳 cost，若代入最佳解 則我們必有
\[
V_T^* \le \bar V_T
\]故可推論 sequence $\{V_T^*\}$ 必定 bounded above。
接著我們證明 optimal cost sequence $\{V_T^*\}$ 為 nondecreasing。
給定量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$，定義 在時間 $T$ 的最佳狀態 sequence 為
\[\left\{ {\hat x\left( {0} \right),\hat x\left( {1} \right),...,\hat x\left( {T} \right)} \right\}
\]現在我們比對 $T$ 時刻的 optimal cost 與 $T-1$ 時刻的 optimal cost 可知
\[V_T^* - \frac{1}{2}\left( {\left| {y\left( {T - 1} \right) - C\hat x\left( {T - 1} \right)} \right|_{{R^{ - 1}}}^2 + |\hat x\left( T \right) - A\hat x\left( k \right)|_{{R^{ - 1}}}^2} \right) \ge V_{T - 1}^*\]此說明了 sequence $\{V_T^*\}$ 為 nondecreasing。
故 由 optimal cost sequence $\{V_T^*\}$ bounded above 與 $\{V_T^*\}$ 為 nondecreasing，我們可推論當 $T \to \infty$ $\{V_T^*\}$ 必定收斂。 $\square$
注意到 optimal estimator cost $V_T^*$ 的收斂性與 系統可觀測性無關，但若我們要求我們的估計狀態 $\hat x \to x$ 則系統的觀測性將扮演重要腳色，我們將此結果記做以下定理：
==========================
Theorem: Estimator Convergence
考慮控制系統 $(A,C)$ 可觀測 且 $Q,R >0$ 為正定矩陣 且 給定一組無雜訊量測輸出
\[
{\bf y}(T) := \{Cx(0), CAx(0), ..., CA^T x(0)\}
\]則 最佳狀態估測 收斂到原本系統狀態；亦即
\[
\hat x(T) \to x(T) \;\; \text{ as $T \to \infty$}
\]==========================
\underbrace {V_{T + n - 1}^* - V_{T - 1}^*}_{ \to 0} \ge \frac{1}{2}\left[ \begin{array}{l}
\sum\limits_{j = - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} \\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right]\\
\Rightarrow \sum\limits_{j = - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2} \to 0\\
\Rightarrow \left\{ \begin{array}{l}
\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j = - 1,...,n - 2\\
y\left( {j + T} \right) - C\hat x\left( {j + T} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j = 0,...,n - 1
\end{array} \right. \ \ \ \ \ \ (**)
\end{array}\]令 $\hat w_T(j) := \hat x(T+j+1|T+n-1) - A \hat x(T+j|T+n-1)$ 並且透過系統方程 $x(k+1) = Ax(k) + w(k)$ 我們有
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I\\
A\\
\vdots \\
{{A^{n - 1}}}
\end{array}} \right]\hat x\left( {T|T + n - {\rm{1}}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
I&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{{A^{n - 2}}}&{{A^{n - 3}}}& \cdots &I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]\ \ \ \ (*)
\end{array} \]且由於我們的量測輸出 滿足
\[\left\{ \begin{array}{l}
y\left( T \right) = Cx\left( T \right)\\
y\left( {T + 1} \right) = CAx\left( T \right)\\
\vdots \\
y\left( {T + n - 1} \right) = C{A^{n - 1}}x\left( T \right)
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] = Ox\left( T \right)\]其中 $O$ 為 observability matrix。現在用上式減去 同乘 $C$ 矩陣 後的 $(*)$ 可得
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]
\end{array}\]現在用 $(**)$ 可知
\[\begin{array}{l}
\underbrace {\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right]}_{ \to 0} = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \underbrace {\left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]}_{ \to 0}\\
\Rightarrow O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right] \to 0
\end{array}\]又因為 observability matrix $O$ 有 linear independent columns，故我們可推知
\[x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right) \to 0 \;\; \text{as $T \to \infty$}\]亦即
\[\hat x\left( {T|T + n - {\rm{1}}} \right) \to x\left( T \right)
\]再者如果我們觀察 $(*)$ 可以發現因為 $\hat w_T (j) \to 0$ 當 $T \to \infty$，故
\[\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right) \;\;\;\; \text{ as $T \to \infty$}\]又因為 $A^{n-1} x(T) = x(T + n -1)$ ，我們有
\[\begin{array}{l}
x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}x\left( T \right) - {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right)\\
\Rightarrow x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\underbrace {\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]}_{ \to 0}\\
\Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)\\
\Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1} \right)\\
\Rightarrow x\left( j \right) \to \hat x\left( j \right),\begin{array}{*{20}{c}}
{}
\end{array}\forall j \to \infty
\end{array}\]至此得證
\[\begin{array}{l}
x(k + 1) = Ax(k) + w(k)\\
y(k) = Cx\left( k \right) + v(k)
\end{array}\]其中 $x \in \mathbb{R}^n, u \in \mathbb{R}^m, y \in \mathbb{R}^p$ 且 $x(0) = x_0$；
基本狀態估計問題：給定初始估計誤差 ( i.e., $\bar x(0) \neq x_0$) 且 不考慮 雜訊 無外部干擾的情況 ($w(k) = v(k)=0$)，我們想問 $\hat x(k) \to x(k)$ as $k \to \infty$ ?
=================
Theorem: (Convergence of Estimator Cost)
給定無 noise 量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$ 則最佳估測器的 cost $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 。
=================
Proof:
由於
\[V_T^{} = \frac{1}{2}\left( \begin{array}{l}
\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{k = 0}^{T - 1} {|\hat x\left( {k + 1} \right) - A\hat x\left( k \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{k = 0}^T {\left| {y\left( k \right) - C\hat x\left( k \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right)\]故我們分三個步驟證明 $V_T^*({\bf y}(T)) $ 在 $T \to \infty$ 時收斂 ：
首先證明 sequence $\{V_T^*\}$ 有上界(bounded above)；接著我們證明 $\{V_T^*\}$ 為 nondecreasing。則由步驟一與步驟二可推論 $\{V_T^*\}$ 必定收斂。
現在我們開始證明 $V_T^*({\bf y}(T)) $ 有界：
由於我們的目標為 $\min_{\hat{ {\bf x}}(k)} V_k$ 但我們並知道 怎樣的 $\hat x$ 可以幫助我們達成此目標，故我們先暫取 $\hat x(0) := x_0$ (並不一定為最佳解!) 則我們有
\[\begin{array}{l}
\left\{ \begin{array}{l}
\hat x(1) = A\hat x(0)\\
\hat x(2) = A\hat x(1) = {A^2}\hat x(0)\\
\vdots
\end{array} \right.\\
\Rightarrow y\left( k \right) = C{A^k}{x_0}
\end{array}\]故若將此結果帶入我們的 cost 可得
\[\bar{V}_T^{} = \frac{1}{2}\left| {\hat x\left( 0 \right) - \bar x\left( 0 \right)} \right|_{{{\left( {{P^ - }\left( 0 \right)} \right)}^{ - 1}}}^2 < \infty
\]但注意到此並非最佳 cost，若代入最佳解 則我們必有
\[
V_T^* \le \bar V_T
\]故可推論 sequence $\{V_T^*\}$ 必定 bounded above。
接著我們證明 optimal cost sequence $\{V_T^*\}$ 為 nondecreasing。
給定量測輸出 ${\bf y}( T) = \{Cx(0), CAx(0),..., CA^T x(0)\}$，定義 在時間 $T$ 的最佳狀態 sequence 為
\[\left\{ {\hat x\left( {0} \right),\hat x\left( {1} \right),...,\hat x\left( {T} \right)} \right\}
\]現在我們比對 $T$ 時刻的 optimal cost 與 $T-1$ 時刻的 optimal cost 可知
\[V_T^* - \frac{1}{2}\left( {\left| {y\left( {T - 1} \right) - C\hat x\left( {T - 1} \right)} \right|_{{R^{ - 1}}}^2 + |\hat x\left( T \right) - A\hat x\left( k \right)|_{{R^{ - 1}}}^2} \right) \ge V_{T - 1}^*\]此說明了 sequence $\{V_T^*\}$ 為 nondecreasing。
故 由 optimal cost sequence $\{V_T^*\}$ bounded above 與 $\{V_T^*\}$ 為 nondecreasing，我們可推論當 $T \to \infty$ $\{V_T^*\}$ 必定收斂。 $\square$
注意到 optimal estimator cost $V_T^*$ 的收斂性與 系統可觀測性無關，但若我們要求我們的估計狀態 $\hat x \to x$ 則系統的觀測性將扮演重要腳色，我們將此結果記做以下定理：
==========================
Theorem: Estimator Convergence
考慮控制系統 $(A,C)$ 可觀測 且 $Q,R >0$ 為正定矩陣 且 給定一組無雜訊量測輸出
\[
{\bf y}(T) := \{Cx(0), CAx(0), ..., CA^T x(0)\}
\]則 最佳狀態估測 收斂到原本系統狀態；亦即
\[
\hat x(T) \to x(T) \;\; \text{ as $T \to \infty$}
\]==========================
Proof:
利用 在時刻 $T+n-1$ 的最佳解作為在時刻 $T-1$的 decision variables ，則前述 Theorem 告訴我們可寫
\[\small V_{T + n - 1}^* - \frac{1}{2}\left[ {\sum\limits_{k = T - 1}^{T + n - 2} {|\hat x\left( {k + 1} \right) - A\hat x\left( k \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{k = T}^{T + n - 1} {\left| {y\left( k \right) - C\hat x\left( k \right)} \right|_{{R^{ - 1}}}^2} } \right] \ge V_{T - 1}^*\]做變數變換 令 $j = k - T$ 我們可得
\[ \small V_{T + n - 1}^* - \frac{1}{2}\left[ {\sum\limits_{j = - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2} } \right] \ge V_{T - 1}^*\]注意到當 $T \to \infty$時， $\{V_T^*\}$ 收斂 且 $Q^{-1}, R^{-1} >0$ 故上式
\[\begin{array}{l}\underbrace {V_{T + n - 1}^* - V_{T - 1}^*}_{ \to 0} \ge \frac{1}{2}\left[ \begin{array}{l}
\sum\limits_{j = - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} \\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2}
\end{array} \right]\\
\Rightarrow \sum\limits_{j = - 1}^{n - 2} {|\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right)|_{{Q^{ - 1}}}^2} + \sum\limits_{j = 0}^{n - 1} {\left| {y\left( {j + T} \right) - C\hat x\left( {j + T} \right)} \right|_{{R^{ - 1}}}^2} \to 0\\
\Rightarrow \left\{ \begin{array}{l}
\hat x\left( {T + j + 1} \right) - A\hat x\left( {T + j} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j = - 1,...,n - 2\\
y\left( {j + T} \right) - C\hat x\left( {j + T} \right) \to 0,\begin{array}{*{20}{c}}
{}
\end{array}\forall j = 0,...,n - 1
\end{array} \right. \ \ \ \ \ \ (**)
\end{array}\]令 $\hat w_T(j) := \hat x(T+j+1|T+n-1) - A \hat x(T+j|T+n-1)$ 並且透過系統方程 $x(k+1) = Ax(k) + w(k)$ 我們有
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I\\
A\\
\vdots \\
{{A^{n - 1}}}
\end{array}} \right]\hat x\left( {T|T + n - {\rm{1}}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
I&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{{A^{n - 2}}}&{{A^{n - 3}}}& \cdots &I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]\ \ \ \ (*)
\end{array} \]且由於我們的量測輸出 滿足
\[\left\{ \begin{array}{l}
y\left( T \right) = Cx\left( T \right)\\
y\left( {T + 1} \right) = CAx\left( T \right)\\
\vdots \\
y\left( {T + n - 1} \right) = C{A^{n - 1}}x\left( T \right)
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] = Ox\left( T \right)\]其中 $O$ 為 observability matrix。現在用上式減去 同乘 $C$ 矩陣 後的 $(*)$ 可得
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right] = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]
\end{array}\]現在用 $(**)$ 可知
\[\begin{array}{l}
\underbrace {\left[ {\begin{array}{*{20}{c}}
{y\left( T \right)}\\
{y\left( {T + 1} \right)}\\
\vdots \\
{y\left( {T + n - 1} \right)}
\end{array}} \right] - C\left[ {\begin{array}{*{20}{c}}
{\hat x\left( {T|T + n - {\rm{1}}} \right)}\\
{\hat x\left( {T{\rm{ + 1}}|T + n - {\rm{1}}} \right)}\\
\vdots \\
{\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)}
\end{array}} \right]}_{ \to 0} = O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} + \underbrace {\left[ {\begin{array}{*{20}{c}}
0&{}&{}&{}\\
C&0&{}&{}\\
\vdots & \vdots & \ddots &{}\\
{C{A^{n - 2}}}&{C{A^{n - 3}}}& \cdots &C
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat w}_T}\left( 0 \right)}\\
{{{\hat w}_T}\left( 1 \right)}\\
\vdots \\
{{{\hat w}_T}\left( {n - 2} \right)}
\end{array}} \right]}_{ \to 0}\\
\Rightarrow O\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right] \to 0
\end{array}\]又因為 observability matrix $O$ 有 linear independent columns，故我們可推知
\[x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right) \to 0 \;\; \text{as $T \to \infty$}\]亦即
\[\hat x\left( {T|T + n - {\rm{1}}} \right) \to x\left( T \right)
\]再者如果我們觀察 $(*)$ 可以發現因為 $\hat w_T (j) \to 0$ 當 $T \to \infty$，故
\[\hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right) \;\;\;\; \text{ as $T \to \infty$}\]又因為 $A^{n-1} x(T) = x(T + n -1)$ ，我們有
\[\begin{array}{l}
x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}x\left( T \right) - {A^{n - 1}}\hat x\left( {T|T + n - {\rm{1}}} \right)\\
\Rightarrow x\left( {T + n - 1} \right) - \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right) \to {A^{n - 1}}\underbrace {\left[ {x\left( T \right) - \hat x\left( {T|T + n - {\rm{1}}} \right)} \right]}_{ \to 0}\\
\Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1|T + n - {\rm{1}}} \right)\\
\Rightarrow x\left( {T + n - 1} \right) \to \hat x\left( {T + n - 1} \right)\\
\Rightarrow x\left( j \right) \to \hat x\left( j \right),\begin{array}{*{20}{c}}
{}
\end{array}\forall j \to \infty
\end{array}\]至此得證
2/19/2015
[控制理論] 離散線性系統的 追蹤 與 調節 問題
Setpoint Tracking
考慮離散線性系統
\[\left\{ \begin{array}{l}
x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\
y\left( k \right) = Cx\left( k \right)
\end{array} \right.\] $x \in \mathbb{R}^n, y \in \mathbb{R}^p, u \in \mathbb{R}^m$；一般而言，控制系統中常見的 追蹤(setpoint tracking)問題 (或稱 servo problem) 如下：
給定 setpoint $y_{sp}$ (e.g., 步階訊號 or 常數)，我們希望在系統達到穩態(steady-state) 之後，系統的輸出 $ y = y_{sp} $。
若我們考慮給定的 setpoint $y_{sp} = 0$ 則稱此類問題為 regulation problem。
回憶在最佳控制理論中的 LQR 方法給予我們對於線性系統的 regulation 問題提供一組最佳解，故我們的想問是否能將此法應用在 一般的追蹤問題?
答案是肯定的，僅需引入 deviation variable 做基本 座標轉換 即可。
現在我們令 $y_{sp}$ 為 output setpoint，定義系統穩態時候的狀態 與 控制力為 $(x_s, u_s)$ 則線性系統方程可改寫為
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_s} = C{x_s}
\end{array} \right.\]注意到對於穩態時 我們希望 $y_s = y_{sp}$ 故 我們有 $y_{sp} = C x_s$ 現在我們將上式改寫為矩陣形式
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]若上述矩陣方程有解 (亦即可以解出 $x_s, u_s$)，且系統控制力無拘束，則我們可以定義 deviation variables 如下
\[\left\{ \begin{array}{l}
\tilde x\left( k \right): = x\left( k \right) - {x_s}\\
\tilde u\left( k \right): = u\left( k \right) - {u_s}
\end{array} \right.\]現在觀察
\[\begin{array}{l}
\tilde x\left( {k + 1} \right) = x\left( {k + 1} \right) - {x_s}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = Ax\left( k \right) + Bu\left( k \right) - \left( {A{x_s} + B{u_s}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = A\left( {x\left( k \right) - {x_s}} \right) + B\left( {u\left( k \right) - {u_s}} \right)\\
\Rightarrow \tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]此表示我們的 deviation variable 仍然滿足原本給定的線性系統。且此時若考慮使用 deviation variable $\tilde{x}, \tilde{u}$ 則 我們的 setpoint tracking problem 被改寫為 regulation problem ；也就是說我們要找到一組 $\tilde{u}(k)$ 使得 $\tilde x(k) \to 0$ (此等價為 $x(k) \to x_s$，故 $C x(k) \to C x_s = y_{sp}$)。解完此 regulation problem 之後，真實的控制力 $u(k) = \tilde u(k) + u_s$。
Comment:
注意到上述論述建立在 $(n + p) \times (n + m)$ 的矩陣方程
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]有解。但何時才有解/與此解是否唯一的問題 我們並未解決 !! 以下我們將對此點進行討論。
對任意 setpoint $y_{sp}$， 上述矩陣方程 $(*)$ 有解 的充分條件為：$(*)$ 要有 linearly independent rows；亦即 $p \le m$ 也就是說我們需要 控制力的數目 $m$ 與 量測輸出的數目 $p$ 至少要相等)。但在實際情況上卻是常常相反，我們可能會獲得非常多的量測輸出數目 (因為裝了很多 sensor)，但實際可以調控的變數卻少於 sensor 數目。為了要解決此問題，我們可選定某矩陣 $H$ 並引入新的 控制變數 $r \in \mathbb{R}^{n_c}$ 為 量測輸出的線性組合；亦即
\[
r := Hy
\]如此一來，若 $p > m$ 情況發生時，我們可選一部份的 輸出 $n_c \le m$ 作為控制變數 並且 對此引入的控制變數 $r$ 也給定所需的 setpoint $r_{sp}$。
另外若 $m > p$ (控制力數目 大於 量測輸出的數目) 時，則對某些 $H$ 與 $r_{sp}$，矩陣方程 $(*)$ 有解 但此時唯一性並不被保證，故我們會希望有唯一解，此時需要對 穩態控制力 $u_s$ 也給定 setpoint $u_{sp}$。
總和以上所述，我們可以建構以下 Steady-State Target Problem
========================
Steady-State Target Problem
考慮最佳化問題
\[
\min_{x_s, u_s} \frac{1}{2} (|u_s - u_{sp}|_{R_s}^2 + |y_s - y_{sp}|_{Q_s}^2)
\]subject to
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right]
\]與 $E u_s \le e$ 與 $F C x_s \le f$。其中 $R_s$ 為 正定矩陣。且對 控制變數的 setpoints $r_{sp}$，上述 target problem 為 feasible 。
========================
有了 steady-state 的解之後我們可以便可以求解 當初我們引入 deviation variable 的 regulation problem:
========================
Dynamic Regulation Problem
考慮以下 cost function
\[\begin{array}{l}
V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right) = \frac{1}{2}\sum\limits_{k = 0}^{N - 1} {\left| {\tilde x\left( k \right)} \right|_Q^2 + \left| {\tilde u\left( k \right)} \right|_R^2} \\
s.t.\\
\tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]其中 $\tilde x(0) = \hat x(k) - x_s$ (此初始條件表示 透過 steady-state $x_s$ 平移估計狀態 $\hat x$ 而得。)；且所求得的 regulator 將會求解以下的 regulation problem
\[\begin{array}{l}
\mathop {\min }\limits_{{\bf{\tilde u}}} V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right)\\
s.t.\\
E\tilde u \le e - E{u_s}\\
FC\tilde x \le f - FC{x_s}
\end{array}\]===================
上述 regulation problem 的 optimal cost 為 $V^*(\tilde x(0))$ 且對應的控制力為 $\tilde {u}^0 (\tilde{x}(0))$
ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design"
考慮離散線性系統
\[\left\{ \begin{array}{l}
x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\
y\left( k \right) = Cx\left( k \right)
\end{array} \right.\] $x \in \mathbb{R}^n, y \in \mathbb{R}^p, u \in \mathbb{R}^m$；一般而言，控制系統中常見的 追蹤(setpoint tracking)問題 (或稱 servo problem) 如下：
給定 setpoint $y_{sp}$ (e.g., 步階訊號 or 常數)，我們希望在系統達到穩態(steady-state) 之後，系統的輸出 $ y = y_{sp} $。
若我們考慮給定的 setpoint $y_{sp} = 0$ 則稱此類問題為 regulation problem。
回憶在最佳控制理論中的 LQR 方法給予我們對於線性系統的 regulation 問題提供一組最佳解，故我們的想問是否能將此法應用在 一般的追蹤問題?
答案是肯定的，僅需引入 deviation variable 做基本 座標轉換 即可。
現在我們令 $y_{sp}$ 為 output setpoint，定義系統穩態時候的狀態 與 控制力為 $(x_s, u_s)$ 則線性系統方程可改寫為
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_s} = C{x_s}
\end{array} \right.\]注意到對於穩態時 我們希望 $y_s = y_{sp}$ 故 我們有 $y_{sp} = C x_s$ 現在我們將上式改寫為矩陣形式
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]若上述矩陣方程有解 (亦即可以解出 $x_s, u_s$)，且系統控制力無拘束，則我們可以定義 deviation variables 如下
\[\left\{ \begin{array}{l}
\tilde x\left( k \right): = x\left( k \right) - {x_s}\\
\tilde u\left( k \right): = u\left( k \right) - {u_s}
\end{array} \right.\]現在觀察
\[\begin{array}{l}
\tilde x\left( {k + 1} \right) = x\left( {k + 1} \right) - {x_s}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = Ax\left( k \right) + Bu\left( k \right) - \left( {A{x_s} + B{u_s}} \right)\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = A\left( {x\left( k \right) - {x_s}} \right) + B\left( {u\left( k \right) - {u_s}} \right)\\
\Rightarrow \tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]此表示我們的 deviation variable 仍然滿足原本給定的線性系統。且此時若考慮使用 deviation variable $\tilde{x}, \tilde{u}$ 則 我們的 setpoint tracking problem 被改寫為 regulation problem ；也就是說我們要找到一組 $\tilde{u}(k)$ 使得 $\tilde x(k) \to 0$ (此等價為 $x(k) \to x_s$，故 $C x(k) \to C x_s = y_{sp}$)。解完此 regulation problem 之後，真實的控制力 $u(k) = \tilde u(k) + u_s$。
Comment:
注意到上述論述建立在 $(n + p) \times (n + m)$ 的矩陣方程
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right] \ \ \ \ (*)
\]有解。但何時才有解/與此解是否唯一的問題 我們並未解決 !! 以下我們將對此點進行討論。
對任意 setpoint $y_{sp}$， 上述矩陣方程 $(*)$ 有解 的充分條件為：$(*)$ 要有 linearly independent rows；亦即 $p \le m$ 也就是說我們需要 控制力的數目 $m$ 與 量測輸出的數目 $p$ 至少要相等)。但在實際情況上卻是常常相反，我們可能會獲得非常多的量測輸出數目 (因為裝了很多 sensor)，但實際可以調控的變數卻少於 sensor 數目。為了要解決此問題，我們可選定某矩陣 $H$ 並引入新的 控制變數 $r \in \mathbb{R}^{n_c}$ 為 量測輸出的線性組合；亦即
\[
r := Hy
\]如此一來，若 $p > m$ 情況發生時，我們可選一部份的 輸出 $n_c \le m$ 作為控制變數 並且 對此引入的控制變數 $r$ 也給定所需的 setpoint $r_{sp}$。
另外若 $m > p$ (控制力數目 大於 量測輸出的數目) 時，則對某些 $H$ 與 $r_{sp}$，矩陣方程 $(*)$ 有解 但此時唯一性並不被保證，故我們會希望有唯一解，此時需要對 穩態控制力 $u_s$ 也給定 setpoint $u_{sp}$。
總和以上所述，我們可以建構以下 Steady-State Target Problem
========================
Steady-State Target Problem
考慮最佳化問題
\[
\min_{x_s, u_s} \frac{1}{2} (|u_s - u_{sp}|_{R_s}^2 + |y_s - y_{sp}|_{Q_s}^2)
\]subject to
\[\left\{ \begin{array}{l}
{x_s} = A{x_s} + B{u_s}\\
{y_{sp}} = C{x_s}
\end{array} \right. \Rightarrow \left[ {\begin{array}{*{20}{c}}
{I - A}&{ - B}\\
C&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_s}}\\
{{u_s}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
0\\
{{y_{sp}}}
\end{array}} \right]
\]與 $E u_s \le e$ 與 $F C x_s \le f$。其中 $R_s$ 為 正定矩陣。且對 控制變數的 setpoints $r_{sp}$，上述 target problem 為 feasible 。
========================
有了 steady-state 的解之後我們可以便可以求解 當初我們引入 deviation variable 的 regulation problem:
========================
Dynamic Regulation Problem
考慮以下 cost function
\[\begin{array}{l}
V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right) = \frac{1}{2}\sum\limits_{k = 0}^{N - 1} {\left| {\tilde x\left( k \right)} \right|_Q^2 + \left| {\tilde u\left( k \right)} \right|_R^2} \\
s.t.\\
\tilde x\left( {k + 1} \right) = A\tilde x\left( k \right) + B\tilde u\left( k \right)
\end{array}\]其中 $\tilde x(0) = \hat x(k) - x_s$ (此初始條件表示 透過 steady-state $x_s$ 平移估計狀態 $\hat x$ 而得。)；且所求得的 regulator 將會求解以下的 regulation problem
\[\begin{array}{l}
\mathop {\min }\limits_{{\bf{\tilde u}}} V\left( {\tilde x\left( 0 \right),{\bf{\tilde u}}} \right)\\
s.t.\\
E\tilde u \le e - E{u_s}\\
FC\tilde x \le f - FC{x_s}
\end{array}\]===================
上述 regulation problem 的 optimal cost 為 $V^*(\tilde x(0))$ 且對應的控制力為 $\tilde {u}^0 (\tilde{x}(0))$
2/15/2015
[線性系統] Hautus Lemma 與 控制性/可穩定性質
以下我們介紹 線性系統理論中關於 控制性的 重要結果
==========================
Lemma: Hautus Lemma for Controllability
我們稱一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 controllable 或稱 (A, B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in \mathbb{C}
\]其中 $\mathbb{C}$ 表示 任意 complex number 所形成的集合。
==========================
Comments:
1. 上述 Lemma 等價為 矩陣 $[\lambda I - A\;\; B]$ 有 $n$ 個 independent rows (full row rank) 。
2. 也許讀者會好奇為何需要如此多 控制性 檢驗工具? 為什麼不使用 可控性矩陣(controllability matrix) 檢驗法就好? 事實上此 Hautus Lemma 主要是功用是大幅簡化 理論證明，但一般實際 檢驗 系統的 可控性 仍多仰賴 可控性矩陣(controllability matrix) 檢驗法。
3. 注意到上述條件需要檢驗 "任意" complex number $\lambda$ 此條件明顯過於嚴苛。不過我們可以注意到 若 $\lambda $ 並非為 $A$ 矩陣的特徵值，則 $rank[\lambda I - A\;\; B] =n$ 必然成立，故我們只需檢驗 $\lambda$ 為 $A$ 矩陣的特徵值部分即可。故
==========================
Lemma (Modified): Hautus Lemma for Controllability
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\]為 controllable 或稱 (A,B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in eig(A)
\]==========================
以下我們將討論拓展到 stabilizability。若現在我們將原本系統 $x(k+1) = A x(k) + Bu(k)$ 透過similar transformation (關於相似轉換細節請參閱 [線性系統] 控制性矩陣 與 非奇異轉換 (Controllability matrix & Non-singular transformation)) 轉成以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。上式又稱為 controllability canonical form。注意到若我們觀察此系統的 contorllability matrix
\[C = \left[ {\begin{array}{*{20}{c}}
{{B_1}}&{{A_{11}}{B_1}}&{{A_{11}}^2{B_1}}& \cdots &{{A_{11}}^n{B_1}}\\
0&0&0& \cdots &0
\end{array}} \right]\]由於 $(A_{11}, B_1)$ 為 controllable 故 controllability matrix $C$ 上排具有 $n_1$ 個 independent row。但是 $C$ 的 下排 $n_2$ rows 均為 $0$ ，故必定不滿足 full row rank test，此系統 uncontrollable。
Exercise: 試證 uncontrollable mode 為 $x_2(k)$。
儘管有 uncontrollable mode ，但我們可以退而求其次，若此 uncontrollable mode 為 stable，則我們仍可控制此系統，此類系統稱作 stabilizable。
===========
Definition: Stabilizability
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。若 $A_{22}$ 為 stable 則此 partitioned system 稱作 stabilizable。=========
接著我們可以拓展前述討論的 Hautus Lemma 到 Stabilizability 之中。
================
Lemma: (Hautus Lemma for Stabilizability)
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 stabilizable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1
\]================
Proof:
$(\Rightarrow)$ 假設系統 stabilizable，我們要證明 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立。
注意到 系統 stabilizable，故由 stabilizability 定義可知
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 。故我們有
\[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right] \ \ \ \ \ (*)
\] 且 若我們觀察 矩陣 $[\lambda I - A_{11}\;\; B_1]$ 與 $[\lambda I - A_{22}]$ 的 row，可發現若 $|\lambda| \ge 1$ 則這些 rows 互為 independent。故由 Hautus Lemma for controllability 可知 ，對 $|\lambda| \ge 1$而言，$(*)$ 的 rows 為 independent 故 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $
$(\Leftarrow)$ 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立，我們要證明 系統 stabilizable 。亦即 考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 注意 partitioned system 為 controllability canonical form，故 $(A_{11}, B_1)$ 已為 controllable ，故我們僅須證明 $A_{22}$ 為 stable 。
由 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $ 我們可知 對任意 $|\lambda| \ge 1 $，
\[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right]\]具有 independent rows 。此暗示了 上述矩陣下排 $[\lambda I - A_{22} ]$ 亦有 full row rank $\forall |\lambda| \ge 1$ 故 $A_{22}$ 的 eigenvalue 必 $<1$ 故 $A_{22} $ 為 stable。
==========================
Lemma: Hautus Lemma for Controllability
我們稱一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 controllable 或稱 (A, B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in \mathbb{C}
\]其中 $\mathbb{C}$ 表示 任意 complex number 所形成的集合。
==========================
Proof. omitted.
1. 上述 Lemma 等價為 矩陣 $[\lambda I - A\;\; B]$ 有 $n$ 個 independent rows (full row rank) 。
2. 也許讀者會好奇為何需要如此多 控制性 檢驗工具? 為什麼不使用 可控性矩陣(controllability matrix) 檢驗法就好? 事實上此 Hautus Lemma 主要是功用是大幅簡化 理論證明，但一般實際 檢驗 系統的 可控性 仍多仰賴 可控性矩陣(controllability matrix) 檢驗法。
3. 注意到上述條件需要檢驗 "任意" complex number $\lambda$ 此條件明顯過於嚴苛。不過我們可以注意到 若 $\lambda $ 並非為 $A$ 矩陣的特徵值，則 $rank[\lambda I - A\;\; B] =n$ 必然成立，故我們只需檢驗 $\lambda$ 為 $A$ 矩陣的特徵值部分即可。故
==========================
Lemma (Modified): Hautus Lemma for Controllability
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\]為 controllable 或稱 (A,B) controllable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall \lambda \in eig(A)
\]==========================
Proof: omitted.
以下我們將討論拓展到 stabilizability。若現在我們將原本系統 $x(k+1) = A x(k) + Bu(k)$ 透過similar transformation (關於相似轉換細節請參閱 [線性系統] 控制性矩陣 與 非奇異轉換 (Controllability matrix & Non-singular transformation)) 轉成以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。上式又稱為 controllability canonical form。注意到若我們觀察此系統的 contorllability matrix
\[C = \left[ {\begin{array}{*{20}{c}}
{{B_1}}&{{A_{11}}{B_1}}&{{A_{11}}^2{B_1}}& \cdots &{{A_{11}}^n{B_1}}\\
0&0&0& \cdots &0
\end{array}} \right]\]由於 $(A_{11}, B_1)$ 為 controllable 故 controllability matrix $C$ 上排具有 $n_1$ 個 independent row。但是 $C$ 的 下排 $n_2$ rows 均為 $0$ ，故必定不滿足 full row rank test，此系統 uncontrollable。
Exercise: 試證 uncontrollable mode 為 $x_2(k)$。
儘管有 uncontrollable mode ，但我們可以退而求其次，若此 uncontrollable mode 為 stable，則我們仍可控制此系統，此類系統稱作 stabilizable。
===========
Definition: Stabilizability
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable。若 $A_{22}$ 為 stable 則此 partitioned system 稱作 stabilizable。=========
接著我們可以拓展前述討論的 Hautus Lemma 到 Stabilizability 之中。
================
Lemma: (Hautus Lemma for Stabilizability)
一個離散線性系統
\[
x(k+1) = A x(k) + Bu(k)
\] 為 stabilizable 若且唯若
\[
rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1
\]================
Proof:
$(\Rightarrow)$ 假設系統 stabilizable，我們要證明 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立。
注意到 系統 stabilizable，故由 stabilizability 定義可知
考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 。故我們有
\[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right] \ \ \ \ \ (*)
\] 且 若我們觀察 矩陣 $[\lambda I - A_{11}\;\; B_1]$ 與 $[\lambda I - A_{22}]$ 的 row，可發現若 $|\lambda| \ge 1$ 則這些 rows 互為 independent。故由 Hautus Lemma for controllability 可知 ，對 $|\lambda| \ge 1$而言，$(*)$ 的 rows 為 independent 故 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $
$(\Leftarrow)$ 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $成立，我們要證明 系統 stabilizable 。亦即 考慮以下 partitioned system
\[\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( {k + 1} \right)}\\
{{x_2}\left( {k + 1} \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{A_{11}}}&{{A_{12}}}\\
0&{{A_{22}}}
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_1}\left( k \right)}\\
{{x_2}\left( k \right)}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
{{B_1}}\\
0
\end{array}} \right]u\left( k \right)\]其中 $(A_{11}, B_1)$ 為 controllable 且 $A_{22}$ 為 stable 注意 partitioned system 為 controllability canonical form，故 $(A_{11}, B_1)$ 已為 controllable ，故我們僅須證明 $A_{22}$ 為 stable 。
由 假設 $rank[\lambda I - A\;\; B] =n, \;\; \forall |\lambda| \ge 1 $ 我們可知 對任意 $|\lambda| \ge 1 $，
\[rank[\lambda I - A\;\;B] = rank\left[ {\begin{array}{*{20}{c}}
{\lambda I - {A_{11}}}&{ - {A_{12}}}&{{B_1}}\\
0&{\lambda I - {A_{22}}}&0
\end{array}} \right]\]具有 independent rows 。此暗示了 上述矩陣下排 $[\lambda I - A_{22} ]$ 亦有 full row rank $\forall |\lambda| \ge 1$ 故 $A_{22}$ 的 eigenvalue 必 $<1$ 故 $A_{22} $ 為 stable。
2/10/2015
[最佳控制] 離散時間 穩態 LQR 控制問題 (2) - LQR convergence
回憶離散時間穩態LQR問題： 定義 cost function
\[\begin{array}{l}
V\left( {x,{\bf{u}}} \right): = \frac{1}{2}\sum\limits_{k = 0}^\infty {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)} \\
s.t.\\
x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\
x\left( 0 \right) = x
\end{array}\]且 $Q, R >0$ 為正定矩陣。若 $(A,B)$ 為可控制 則 最佳化問題
\[
\min_{\bf u} V(x,{\bf u})
\]解存在 且唯一 (why?)。
現在我們將此解 (控制力序列) 記做 ${\bf u}^*(x)$；其中第一項為 $u^*(x)$。且對應的最佳控制律 $K_{\infty}(\cdot)$ 故我們有
\[
K_\infty(x) = u^*(x) = {\bf u}^* (0;x)
\]
以下定理說明若可控制系統 採用 穩態 LQR 控制，則保證閉迴路穩定度：
=====================
Theorem: Steady State LQR Guarantees the Closed-Loop Stability
對 $(A,B)$ 可控制，穩態LQR問題 且 $Q,R >0$ 保證閉迴路系統
\[
x(k+1) = A x(k) + B K_\infty (x)
\]穩定。
=====================
Proof:
首先證明 穩態 cost function 有界。注意到由於 $(A,B)$ 可控制，由定義可知存在 一組控制力
\[
{\bf u} = \{u(0), u(1),...,u(n-1)\}
\] 使得 控制後的系統能在有限時間 $n$ 步之內，可由任意給定初始狀態 $x(0)$ 移動到 給定的終止狀態 $x(n) = 0$；且 在時間 $ n$ 之後的控制力 $ \{u(n+1), u(n+2),...\} = \{0,0,0,...\}$ 故原本無窮和的 cost 函數變為有限和
\[\begin{array}{l}
V\left( {x,{\bf{u}}} \right): = \frac{1}{2}\sum\limits_{k = 0}^\infty {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)} \\
\Rightarrow V\left( {x,{\bf{u}}} \right) = \frac{1}{2}\sum\limits_{k = 0}^n {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)}
\end{array}\]且由於 成本函數為 strictly convex in $\bf u$ 且 $R>0$ ( $u$ 沒有 vanish) 故 此穩態LQR最佳問題的解為唯一。
現在我們觀察 costs to go 的數列滿足
\[
V_{k+1} = V_k - 1/2 (x(k)'Qx(k) + u(k)' R u(k))
\]其中 $V_k := V^0(x(k))$ 為 在時刻 $k$ 對應 狀態為 $x(k)$ 與 對應最佳控制力 $u(k) = u^0(x(k))$ 的 cost。故此數列 $\{V_k\}$ 為非遞增 且 有下界 (為 $0$)。可推知此數列必定收斂；故
\[
|V_{k+1} - V_k| \to 0\;\;\;\; \text{as $k \to \infty$}
\]因此
\[\left\{ \begin{array}{l}
x\left( k \right)'Qx\left( k \right)\; \to 0\\
u\left( k \right)'{\rm{ }}R{\rm{ }}u\left( k \right) \to 0
\end{array} \right.\]又因為 $Q,R >0$ 故可推知
\[\left\{ \begin{array}{l}
x\left( k \right)\; \to 0\\
u\left( k \right) \to 0
\end{array} \right.\]亦即閉迴路系統狀態收斂 (閉迴路穩定!)。 $\square$
Comments
1. $R>0$ 條件用以保證控制律唯一性。實際上若使用者不關心控制律，可將其設為參數非常小的正定矩陣
2. 上述定裡假設可從 controllability 推廣到 stabilizability。
3. $Q>0$ 的條件可被推廣到 $Q \ge 0$ 與 $(A,Q)$ detectable。
\[\begin{array}{l}
V\left( {x,{\bf{u}}} \right): = \frac{1}{2}\sum\limits_{k = 0}^\infty {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)} \\
s.t.\\
x\left( {k + 1} \right) = Ax\left( k \right) + Bu\left( k \right)\\
x\left( 0 \right) = x
\end{array}\]且 $Q, R >0$ 為正定矩陣。若 $(A,B)$ 為可控制 則 最佳化問題
\[
\min_{\bf u} V(x,{\bf u})
\]解存在 且唯一 (why?)。
現在我們將此解 (控制力序列) 記做 ${\bf u}^*(x)$；其中第一項為 $u^*(x)$。且對應的最佳控制律 $K_{\infty}(\cdot)$ 故我們有
\[
K_\infty(x) = u^*(x) = {\bf u}^* (0;x)
\]
以下定理說明若可控制系統 採用 穩態 LQR 控制，則保證閉迴路穩定度：
=====================
Theorem: Steady State LQR Guarantees the Closed-Loop Stability
對 $(A,B)$ 可控制，穩態LQR問題 且 $Q,R >0$ 保證閉迴路系統
\[
x(k+1) = A x(k) + B K_\infty (x)
\]穩定。
=====================
Comment:
線性系統中，漸進穩定度 等價 漸進收斂 $x(k) \to 0$
線性系統中，漸進穩定度 等價 漸進收斂 $x(k) \to 0$
Proof:
首先證明 穩態 cost function 有界。注意到由於 $(A,B)$ 可控制，由定義可知存在 一組控制力
\[
{\bf u} = \{u(0), u(1),...,u(n-1)\}
\] 使得 控制後的系統能在有限時間 $n$ 步之內，可由任意給定初始狀態 $x(0)$ 移動到 給定的終止狀態 $x(n) = 0$；且 在時間 $ n$ 之後的控制力 $ \{u(n+1), u(n+2),...\} = \{0,0,0,...\}$ 故原本無窮和的 cost 函數變為有限和
\[\begin{array}{l}
V\left( {x,{\bf{u}}} \right): = \frac{1}{2}\sum\limits_{k = 0}^\infty {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)} \\
\Rightarrow V\left( {x,{\bf{u}}} \right) = \frac{1}{2}\sum\limits_{k = 0}^n {{x^T}\left( k \right)Qx\left( k \right) + {u^T}\left( k \right)Ru\left( k \right)}
\end{array}\]且由於 成本函數為 strictly convex in $\bf u$ 且 $R>0$ ( $u$ 沒有 vanish) 故 此穩態LQR最佳問題的解為唯一。
現在我們觀察 costs to go 的數列滿足
\[
V_{k+1} = V_k - 1/2 (x(k)'Qx(k) + u(k)' R u(k))
\]其中 $V_k := V^0(x(k))$ 為 在時刻 $k$ 對應 狀態為 $x(k)$ 與 對應最佳控制力 $u(k) = u^0(x(k))$ 的 cost。故此數列 $\{V_k\}$ 為非遞增 且 有下界 (為 $0$)。可推知此數列必定收斂；故
\[
|V_{k+1} - V_k| \to 0\;\;\;\; \text{as $k \to \infty$}
\]因此
\[\left\{ \begin{array}{l}
x\left( k \right)'Qx\left( k \right)\; \to 0\\
u\left( k \right)'{\rm{ }}R{\rm{ }}u\left( k \right) \to 0
\end{array} \right.\]又因為 $Q,R >0$ 故可推知
\[\left\{ \begin{array}{l}
x\left( k \right)\; \to 0\\
u\left( k \right) \to 0
\end{array} \right.\]亦即閉迴路系統狀態收斂 (閉迴路穩定!)。 $\square$
Comments
1. $R>0$ 條件用以保證控制律唯一性。實際上若使用者不關心控制律，可將其設為參數非常小的正定矩陣
2. 上述定裡假設可從 controllability 推廣到 stabilizability。
3. $Q>0$ 的條件可被推廣到 $Q \ge 0$ 與 $(A,Q)$ detectable。
2/07/2015
[最佳控制] 線性系統的最佳參數估計 - Kalman Filter (1)
此文章基於前篇 [最佳控制] 線性系統的最佳參數估計 - Kalman Filter (0)；強烈建議讀者先參閱前篇再行閱讀此文。此文將利用 計算 probability density 方式來求解 Kalman filter。
考慮離散時間動態系統
\[
x(k+1) = Ax(k) + w(k);\;\; y(k) = Cx(k) + v(k)
\]且 $w \sim N(0,Q)$ 為製程雜訊 process noise 或稱 干擾 process disturbance；$v \sim N(0,R)$ 為量測雜訊 (measurement noise) ； $x(0) \sim N(\bar x(0), Q(0))$；$x \in \mathbb{R}^n, A \in \mathbb{R}^{n \times n}, C \in \mathbb{R}^{p \times n}, y \in \mathbb{R}^p$，Kalman filter 便是在試圖回答：假設狀態未知我們只能拿到量測輸出，則最佳的狀態估計該是如何?
First Step : $k=0$
假設 初始狀態 $x(0)$ 為 mean $\bar x(0)$ 且 convariance matrix $Q(0)$ 的 normal distrbuted 隨機向量，亦即
\[
x(0) \sim N(\bar x(0), Q(0))
\]接著獲得 初始 (受雜訊污染的) 量測輸出 $y(0)$ 滿足下式
\[
y(0) = C x(0) + v(0)
\]其中 $v(0) \sim N(0, R)$ 為 雜訊 (measurement noise)。
我們的目標：獲得 conditional density $p_{x(0)|y(0)}(x(0) | y(0))$ ，則我們的狀態估計 $\hat x$ 即可透過此 conditional density 求得
$$\hat x := \arg \max_x p_{x(0)|y(0)}(x(0) | y(0))$$
Comments:
1. 若未知 $\bar x(0)$ 或者 $Q(0)$ 則我們通常選 $\bar x(0) :=0$ 且 $Q(0)$ 很大 來表示我們對初始狀態所知甚少 (noninformative prior)。
2. 若量測過程中受到大的雜訊，則我們會給予較大的 $R$。若量測過程十分精確沒有太多雜訊汙染我們的輸出 $y$ 則 $R$ 較小。
3. Conditional density 描述了在我們獲得 初始量測輸出 $y(0)$ 之後我們對於 $x(0)$ 的了解。
現在回歸我們的目標，究竟該如何推得 $p_{x(0)|y(0)}(x(0) | y(0))$ ?
首先考慮 $(x(0), y(0))$ 如下
\[\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{v\left( 0 \right)}
\end{array}} \right]
\]假設 $v(0)$ 與 $x(0)$ 彼此互為獨立。注意到 $x(0), v(0)$ 為 joint normal 亦即
\[\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{v\left( 0 \right)}
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&{{R}}
\end{array}} \right]} \right)\]
上述 $[x(0) \;\; y(0)]$ 為 $[x(0) \;\; v(0)]$ 線性轉換 且 故我們可以馬上知道
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&R
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
{C\bar x\left( 0 \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&R
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
I&{{C^T}}\\
0&I
\end{array}} \right]} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
{C\bar x\left( 0 \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&{Q\left( 0 \right){C^T}}\\
{CQ\left( 0 \right)}&{CQ\left( 0 \right){C^T} + R}
\end{array}} \right]} \right)
\end{array}\]現在有了 $x(0), y(0)$ 的 joint density， 注意到此 joint density 為 normal，故若要計算 conditonal density of $x(0)$ given $y(0)$，亦即 $p_{x(0)|y(0)}(x(0)|y(0))$ ；則我們可以使用前述文章討論的 FACT 3 來求得，亦即
\[
p_{x(0)|y(0)}(x(0)|y(0)) = n(x(0), m, P)
\]其中
\[\left\{ {\begin{array}{*{20}{l}}
\begin{array}{l}
m = \bar x\left( 0 \right) + L\left( 0 \right)(y\left( 0 \right) - C\bar x\left( 0 \right))\\
L\left( 0 \right): = Q\left( 0 \right){C^T}\left( {CQ\left( 0 \right){C^T} + R} \right)_{}^{ - 1}
\end{array}\\
{P = Q\left( 0 \right) - Q\left( 0 \right){C^T}\left( {CQ\left( 0 \right){C^T} + R} \right)_{}^{ - 1}CQ\left( 0 \right)}
\end{array}} \right.\]則 optimal state estimation $\hat x$ 及為 具有最大 conditional density 的 $x(0)$；對於 normal distribution 而言，此 $\hat x$ 剛好為 mean；故我們選 $\hat x(0) := m$；且上式中 $P := P(0)$ 表示獲得 量測輸出 $y(0)$ 之後的 variance
Next Step: State Prediction
考慮現在狀態從 $k=0$ 移動到 $k=1$，我們有
\[
x(1) = Ax(0) + w(0)
\]亦可將其寫成線性轉換型式：
\[x\left( 1 \right) = \left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]\]其中 $w(0) \sim N(0, Q)$ 為 系統干擾 (disturbance) 或稱 製程雜訊 (process noise)。若 狀態受到大的外部擾動，則可以想見 會有較大的 $Q$ convaraince matrix。同理若外部干擾較小則有較小的雜訊。
故我們目標：要計算 conditional density $p_{x(1)|y(0)}(x(1),y(0))$
現在我們需要 conditional on joint density $(x(0),w(0))$ given $y(0)$，假設 $w(0)$ 與 $x(0),v(0)$ 彼此互為獨立，則我們可寫
\[\left( {\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]} \right)\]故現在利用 前述文章討論的 FACT 3' 來求得 conditional density，亦即
\[\begin{array}{l}
\left( {\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]} \right)\\
\Rightarrow \left( {x\left( 1 \right)|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left( {x\left( 1 \right)|y\left( 0 \right)} \right)\sim N\left( {A\hat x\left( 0 \right),AP\left( 0 \right){A^T} + Q} \right)
\end{array}\]故 conditional density 仍為 normal
\[
p_{x(1)|y(0)}(x(1)|y(0)) = n(x(1), \hat x^-(1), P^-(1))
\]其中
\[\left\{ \begin{array}{l}
{{\hat x}^ - }\left( 1 \right) = A\hat x\left( 0 \right)\\
{P^ - }\left( 1 \right) = AP\left( 0 \right){A^T} + Q
\end{array} \right.\]接著我們僅需 遞迴重複上述步驟 $k=2,3,4,...$ 即可。以下我們給出總結：
Summary
定義 量測輸出從初始 直到 時間 $k$ 則
\[
{\bf y}(k) := \{y(0),y(1),...,y(k)\}
\]在 時間 $k$ 時，conditional density with data ${\bf y}(k-1)$ 為 normal
\[
p_{x(k)| {\bf y}(k-1)}(x(k) | {\bf y}(k-1)) = n(x(k), \hat x^-(k), P^-(k))
\]上述 mean 與 covariance matrix 有上標 $^-$ 號表示此估計為僅僅透過 量測 ${\bf y}(k-1)$ 的輸出，並未獲得 $k$ 時刻的量測輸出。 (表示用過去 $k-1$ 資料 預測 $k$ 的狀態! ) 注意到在 $k=0$，遞迴起始於 $\hat x^- (0) = \bar x(0)$ 與 $P^-(0) = Q(0)$。
接著我們會獲得 $y(k)$ 滿足
\[\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]
\]上式表示線性轉換。由於 量測雜訊 $v(k)$ 與 $x(k)$ 以及 ${\bf y}(k-1)$ 彼此獨立，故我們有 density of $(x(k), v(k))$
\[\underbrace {\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]}_{ = \left( {\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]|{\bf{y}}\left( {k - 1} \right)} \right)}\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&0\\
0&R
\end{array}} \right]} \right)\]透過線性轉換結果，可得
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&0\\
0&R
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
{C{{\hat x}^ - }\left( k \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&{{P^ - }\left( k \right){C^T}}\\
{C{P^ - }\left( k \right)}&{C{P^ - }\left( k \right){C^T} + R}
\end{array}} \right]} \right)
\end{array}\]注意到 $\{{\bf y}(k-1), y(k)\} \equiv {\bf y}(k)$ 故利用 conditional density 結果可得
\[
p_{x(k)|{\bf y} (k)}(x(k)| {\bf y}(k)) = n(x(k), \hat x(k), P(k))
\]其中
\[\left\{ {\begin{array}{*{20}{l}}
\begin{array}{l}
\hat x\left( k \right) = {{\hat x}^ - }\left( k \right) + L\left( k \right)(y\left( k \right) - C{{\hat x}^ - }\left( k \right))\\
L\left( k \right) = {P^ - }\left( k \right){C^T}\left( {C{P^ - }\left( k \right){C^T} + R} \right)_{}^{ - 1}
\end{array}\\
{P\left( k \right) = {P^ - }\left( k \right) - {P^ - }\left( k \right){C^T}\left( {C{P^ - }\left( k \right){C^T} + R} \right)_{}^{ - 1}C{P^ - }\left( k \right)}
\end{array}} \right.\]現在我們用下列模型來預估 基於 $k$ 時刻量測值 預估 $k+1$ 時刻
\[x\left( {k + 1} \right) = \left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{w\left( k \right)}
\end{array}} \right]\]由於 $w(k)$ 與 $x(k)$ 以及 $\bf y$$(k)$ 彼此獨立，故 joint density of $(x(k),w(k))$ 可寫為
\[\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{w\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( k \right)}&0\\
0&Q
\end{array}} \right]} \right)\]故 conditional density 為
\[{p_{x\left( {k + 1} \right)|{\bf{y}}\left( k \right)}}\left( {x\left( {k + 1} \right)|{\bf{y}}\left( k \right)} \right) = n\left( {x\left( {k + 1} \right),{{\hat x}^ - }\left( {k + 1} \right),{P^ - }\left( {k + 1} \right)} \right)\]其中
\[\left\{ \begin{array}{l}
{{\hat x}^ - }\left( {k + 1} \right) = A\hat x\left( k \right)\\
{P^ - }\left( {k + 1} \right) = AP\left( k \right){A^T} + Q
\end{array} \right.\]
考慮離散時間動態系統
\[
x(k+1) = Ax(k) + w(k);\;\; y(k) = Cx(k) + v(k)
\]且 $w \sim N(0,Q)$ 為製程雜訊 process noise 或稱 干擾 process disturbance；$v \sim N(0,R)$ 為量測雜訊 (measurement noise) ； $x(0) \sim N(\bar x(0), Q(0))$；$x \in \mathbb{R}^n, A \in \mathbb{R}^{n \times n}, C \in \mathbb{R}^{p \times n}, y \in \mathbb{R}^p$，Kalman filter 便是在試圖回答：假設狀態未知我們只能拿到量測輸出，則最佳的狀態估計該是如何?
First Step : $k=0$
假設 初始狀態 $x(0)$ 為 mean $\bar x(0)$ 且 convariance matrix $Q(0)$ 的 normal distrbuted 隨機向量，亦即
\[
x(0) \sim N(\bar x(0), Q(0))
\]接著獲得 初始 (受雜訊污染的) 量測輸出 $y(0)$ 滿足下式
\[
y(0) = C x(0) + v(0)
\]其中 $v(0) \sim N(0, R)$ 為 雜訊 (measurement noise)。
我們的目標：獲得 conditional density $p_{x(0)|y(0)}(x(0) | y(0))$ ，則我們的狀態估計 $\hat x$ 即可透過此 conditional density 求得
$$\hat x := \arg \max_x p_{x(0)|y(0)}(x(0) | y(0))$$
Comments:
1. 若未知 $\bar x(0)$ 或者 $Q(0)$ 則我們通常選 $\bar x(0) :=0$ 且 $Q(0)$ 很大 來表示我們對初始狀態所知甚少 (noninformative prior)。
2. 若量測過程中受到大的雜訊，則我們會給予較大的 $R$。若量測過程十分精確沒有太多雜訊汙染我們的輸出 $y$ 則 $R$ 較小。
3. Conditional density 描述了在我們獲得 初始量測輸出 $y(0)$ 之後我們對於 $x(0)$ 的了解。
現在回歸我們的目標，究竟該如何推得 $p_{x(0)|y(0)}(x(0) | y(0))$ ?
首先考慮 $(x(0), y(0))$ 如下
\[\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{v\left( 0 \right)}
\end{array}} \right]
\]假設 $v(0)$ 與 $x(0)$ 彼此互為獨立。注意到 $x(0), v(0)$ 為 joint normal 亦即
\[\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{v\left( 0 \right)}
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&{{R}}
\end{array}} \right]} \right)\]
上述 $[x(0) \;\; y(0)]$ 為 $[x(0) \;\; v(0)]$ 線性轉換 且 故我們可以馬上知道
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&R
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
{C\bar x\left( 0 \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&0\\
0&R
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
I&{{C^T}}\\
0&I
\end{array}} \right]} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{y\left( 0 \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\bar x\left( 0 \right)}\\
{C\bar x\left( 0 \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{Q\left( 0 \right)}&{Q\left( 0 \right){C^T}}\\
{CQ\left( 0 \right)}&{CQ\left( 0 \right){C^T} + R}
\end{array}} \right]} \right)
\end{array}\]現在有了 $x(0), y(0)$ 的 joint density， 注意到此 joint density 為 normal，故若要計算 conditonal density of $x(0)$ given $y(0)$，亦即 $p_{x(0)|y(0)}(x(0)|y(0))$ ；則我們可以使用前述文章討論的 FACT 3 來求得，亦即
\[
p_{x(0)|y(0)}(x(0)|y(0)) = n(x(0), m, P)
\]其中
\[\left\{ {\begin{array}{*{20}{l}}
\begin{array}{l}
m = \bar x\left( 0 \right) + L\left( 0 \right)(y\left( 0 \right) - C\bar x\left( 0 \right))\\
L\left( 0 \right): = Q\left( 0 \right){C^T}\left( {CQ\left( 0 \right){C^T} + R} \right)_{}^{ - 1}
\end{array}\\
{P = Q\left( 0 \right) - Q\left( 0 \right){C^T}\left( {CQ\left( 0 \right){C^T} + R} \right)_{}^{ - 1}CQ\left( 0 \right)}
\end{array}} \right.\]則 optimal state estimation $\hat x$ 及為 具有最大 conditional density 的 $x(0)$；對於 normal distribution 而言，此 $\hat x$ 剛好為 mean；故我們選 $\hat x(0) := m$；且上式中 $P := P(0)$ 表示獲得 量測輸出 $y(0)$ 之後的 variance
Next Step: State Prediction
考慮現在狀態從 $k=0$ 移動到 $k=1$，我們有
\[
x(1) = Ax(0) + w(0)
\]亦可將其寫成線性轉換型式：
\[x\left( 1 \right) = \left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]\]其中 $w(0) \sim N(0, Q)$ 為 系統干擾 (disturbance) 或稱 製程雜訊 (process noise)。若 狀態受到大的外部擾動，則可以想見 會有較大的 $Q$ convaraince matrix。同理若外部干擾較小則有較小的雜訊。
故我們目標：要計算 conditional density $p_{x(1)|y(0)}(x(1),y(0))$
現在我們需要 conditional on joint density $(x(0),w(0))$ given $y(0)$，假設 $w(0)$ 與 $x(0),v(0)$ 彼此互為獨立，則我們可寫
\[\left( {\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]} \right)\]故現在利用 前述文章討論的 FACT 3' 來求得 conditional density，亦即
\[\begin{array}{l}
\left( {\left[ {\begin{array}{*{20}{c}}
{x\left( 0 \right)}\\
{w\left( 0 \right)}
\end{array}} \right]|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]} \right)\\
\Rightarrow \left( {x\left( 1 \right)|y\left( 0 \right)} \right)\sim N\left( {\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{\hat x\left( 0 \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{P\left( 0 \right)}&0\\
0&Q
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left( {x\left( 1 \right)|y\left( 0 \right)} \right)\sim N\left( {A\hat x\left( 0 \right),AP\left( 0 \right){A^T} + Q} \right)
\end{array}\]故 conditional density 仍為 normal
\[
p_{x(1)|y(0)}(x(1)|y(0)) = n(x(1), \hat x^-(1), P^-(1))
\]其中
\[\left\{ \begin{array}{l}
{{\hat x}^ - }\left( 1 \right) = A\hat x\left( 0 \right)\\
{P^ - }\left( 1 \right) = AP\left( 0 \right){A^T} + Q
\end{array} \right.\]接著我們僅需 遞迴重複上述步驟 $k=2,3,4,...$ 即可。以下我們給出總結：
Summary
定義 量測輸出從初始 直到 時間 $k$ 則
\[
{\bf y}(k) := \{y(0),y(1),...,y(k)\}
\]在 時間 $k$ 時，conditional density with data ${\bf y}(k-1)$ 為 normal
\[
p_{x(k)| {\bf y}(k-1)}(x(k) | {\bf y}(k-1)) = n(x(k), \hat x^-(k), P^-(k))
\]上述 mean 與 covariance matrix 有上標 $^-$ 號表示此估計為僅僅透過 量測 ${\bf y}(k-1)$ 的輸出，並未獲得 $k$ 時刻的量測輸出。 (表示用過去 $k-1$ 資料 預測 $k$ 的狀態! ) 注意到在 $k=0$，遞迴起始於 $\hat x^- (0) = \bar x(0)$ 與 $P^-(0) = Q(0)$。
接著我們會獲得 $y(k)$ 滿足
\[\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]
\]上式表示線性轉換。由於 量測雜訊 $v(k)$ 與 $x(k)$ 以及 ${\bf y}(k-1)$ 彼此獨立，故我們有 density of $(x(k), v(k))$
\[\underbrace {\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]}_{ = \left( {\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{v\left( k \right)}
\end{array}} \right]|{\bf{y}}\left( {k - 1} \right)} \right)}\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&0\\
0&R
\end{array}} \right]} \right)\]透過線性轉換結果，可得
\[\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&0\\
0&R
\end{array}} \right]{{\left[ {\begin{array}{*{20}{c}}
I&0\\
C&I
\end{array}} \right]}^T}} \right)\\
\Rightarrow \left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{y\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{{\hat x}^ - }\left( k \right)}\\
{C{{\hat x}^ - }\left( k \right)}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P^ - }\left( k \right)}&{{P^ - }\left( k \right){C^T}}\\
{C{P^ - }\left( k \right)}&{C{P^ - }\left( k \right){C^T} + R}
\end{array}} \right]} \right)
\end{array}\]注意到 $\{{\bf y}(k-1), y(k)\} \equiv {\bf y}(k)$ 故利用 conditional density 結果可得
\[
p_{x(k)|{\bf y} (k)}(x(k)| {\bf y}(k)) = n(x(k), \hat x(k), P(k))
\]其中
\[\left\{ {\begin{array}{*{20}{l}}
\begin{array}{l}
\hat x\left( k \right) = {{\hat x}^ - }\left( k \right) + L\left( k \right)(y\left( k \right) - C{{\hat x}^ - }\left( k \right))\\
L\left( k \right) = {P^ - }\left( k \right){C^T}\left( {C{P^ - }\left( k \right){C^T} + R} \right)_{}^{ - 1}
\end{array}\\
{P\left( k \right) = {P^ - }\left( k \right) - {P^ - }\left( k \right){C^T}\left( {C{P^ - }\left( k \right){C^T} + R} \right)_{}^{ - 1}C{P^ - }\left( k \right)}
\end{array}} \right.\]現在我們用下列模型來預估 基於 $k$ 時刻量測值 預估 $k+1$ 時刻
\[x\left( {k + 1} \right) = \left[ {\begin{array}{*{20}{c}}
A&I
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{w\left( k \right)}
\end{array}} \right]\]由於 $w(k)$ 與 $x(k)$ 以及 $\bf y$$(k)$ 彼此獨立，故 joint density of $(x(k),w(k))$ 可寫為
\[\left[ {\begin{array}{*{20}{c}}
{x\left( k \right)}\\
{w\left( k \right)}
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{\hat x\left( k \right)}\\
0
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{P\left( k \right)}&0\\
0&Q
\end{array}} \right]} \right)\]故 conditional density 為
\[{p_{x\left( {k + 1} \right)|{\bf{y}}\left( k \right)}}\left( {x\left( {k + 1} \right)|{\bf{y}}\left( k \right)} \right) = n\left( {x\left( {k + 1} \right),{{\hat x}^ - }\left( {k + 1} \right),{P^ - }\left( {k + 1} \right)} \right)\]其中
\[\left\{ \begin{array}{l}
{{\hat x}^ - }\left( {k + 1} \right) = A\hat x\left( k \right)\\
{P^ - }\left( {k + 1} \right) = AP\left( k \right){A^T} + Q
\end{array} \right.\]
[最佳控制] 線性系統的最佳參數估計 Kalman Filter (0)- 預備知識
這次要介紹 Kalman filter 或稱 Optimal Linear State Estmator，以下我們將簡單介紹一些在下一篇文章需要使用的一些結果。
Preliminary
回憶若 $x$ 為 mean $m$ 且 variance $\sigma^2$ 的 Normal 隨機變數 則 我們表示為 $x \sim N(m, \sigma^2)$ 。其 probability density function 可記做
\[
\frac{1}{\sqrt{2 \pi} \sigma} \exp(\frac{-1}{2 \sigma^2}(x-m)^2)
\]現在若我們拓展上述結果到 多個隨機變數 (又稱 random vector 隨機向量) 的情況，令 $x$ 為Normal 隨機向量 記做
\[
x \sim N(m,P); \;\;
p_x(x) := n(x,m,P)
\]上述符號 表示 $x$ 為 normal distributed 且 mean vector $m$ 與 convariance matrix $P$。另外 $n(x,m,P)$ 表示 normal probability density function
\[\;n\left( {x,m,P} \right): = \frac{1}{{{{\left( {2\pi } \right)}^{n/2}}{{\left( {\det P} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right)\]
Comment:
若 $x \in \mathbb{R}^n$ 則 mean vector $m \in \mathbb{R}^n$ 且 convariance matrix $P \in \mathbb{R}^{n \times n}$ 且為 實數 對稱 正定 矩陣 (正定條件用以確保 $P^{-1}$ 存在，使得上述的 probability density function 可以被定義)。若 $P$ 不為對 正定 我們稱為 singular normal distribution 或稱 degenerate normal。
Example
若 $n=2$ 則我們可以繪製 normal density function； 比如說
\[m = \left[ {\begin{array}{*{20}{c}}
0\\
0
\end{array}} \right];\begin{array}{*{20}{c}}
{}&{}
\end{array}{P^{ - 1}} = \left[ {\begin{array}{*{20}{c}}
{3.5}&{2.5}\\
{2.5}&{4.0}
\end{array}} \right]\]則我們可繪製
現在我們看幾個 之後會使用到的基本結果：
================
FACT 1: Joint independent normals
若隨機向量 $x \sim N(m_x,P_x)$ 與 $y \sim N(m_y,P_y)$ 為 normally distributed 且 彼此互為獨立，則其 joint density $p_{x,y}(x,y)$ 如下
\[
p_{x,y}(x,y) = p_{x}(x)p_{y}(y)=n(x,m_x,P_x) \cdot n(y,m_y,P_y)
\]且
\[\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&0\\
0&{{P_y}}
\end{array}} \right]} \right)\]================
================
FACT 2: Linear Transformation of a normal
若 $x \sim N(m, P)$ 且 $y$ 為 $x$ 的線性轉換；亦即對任意矩陣 $A$， $y = Ax$ 則
\[
y \sim N(Am , APA^T)
\]================
================
FACT 3: Conditional of a joint normal
若 $x,y$ 為 jointly normal distributed (no independent assumption)
\[\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)\]則 conditional density of $x$ given $y$ 仍為 Normal 亦即
\[
(x|y) \sim N(m,P)
\]其 probability density function 為 $ p_{x|y} (x|y) = n(x,m,P)$ 其中 conditional mean vector $m$ 與 conditional convariance matrix $P$ 分別為
\[\begin{array}{l}
m = {m_x} + {P_{xy}}P_y^{ - 1}(y - {m_y})\\
P = {P_x} - {P_{xy}}P_y^{ - 1}{P_{yx}}
\end{array}\]================
注意到上式中 conditional mean $m$ 為 random vector (depends on $y$) 。
Proof:
回憶 conditional density of $x$ given y 定義
\[
p_{x|y}(x,y) := \frac{p_{x,y}(x,y)}{p_y(y)}
\]注意到 $(x,y)$ 為 joint normal 故我們有
\[\small \begin{array}{l}
{p_y}\left( y \right): = \frac{1}{{{{\left( {2\pi } \right)}^{{n_y}/2}}{{\left( {\det {P_y}} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)\\
{p_{x,y}}\left( {x,y} \right): = \frac{1}{{{{\left( {2\pi } \right)}^{\left( {{n_x} + {n_y}} \right)/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]} \right)
\end{array}\]亦即
\[\begin{array}{l}
{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \frac{{\frac{1}{{{{\left( {2\pi } \right)}^{\left( {{n_x} + {n_y}} \right)/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]} \right)}}{{\frac{1}{{{{\left( {2\pi } \right)}^{{n_y}/2}}{{\left( {\det {P_y}} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right] - {{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}
\end{array}\]注意到若我們取 $P:= P_x - P_{xy}P_y^{-1}P_{yx}$ 則 利用 Matrix inversion Lemma 可知
\[{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]^{ - 1}} = \left[ {\begin{array}{*{20}{c}}
{{P^{ - 1}}}&{ - {P^{ - 1}}{P_{xy}}P_y^{ - 1}}\\
{ - P_y^{ - 1}{P_{yx}}{P^{ - 1}}}&{P_y^{ - 1} + P_y^{ - 1}{P_{yx}}{P^{ - 1}}{P_{xy}}P_y^{ - 1}}
\end{array}} \right]\]將此結果帶回我們可得
\[\begin{array}{l}
{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( \begin{array}{l}
{\left( {x - {m_x}} \right)^T}{P^{ - 1}}\left( {x - {m_x}} \right) - 2{\left( {y - {m_y}} \right)^T}P_y^{ - 1}{P_{yx}}{P^{ - 1}}\left( {x - {m_x}} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + {\left( {y - {m_y}} \right)^T}\left( {P_y^{ - 1}{P_{yx}}{P^{ - 1}}{P_{xy}}P_y^{ - 1}} \right)\left( {y - {m_y}} \right)
\end{array} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {{{\left( {x - {m_x}} \right)}^T} - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]}^T}{P^{ - 1}}\left[ {{{\left( {x - {m_x}} \right)}^T} - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {\left( {x - {m_x}} \right) - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]}^T}{P^{ - 1}}\left[ {\left( {x - {m_x}} \right) - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}
\end{array}\]接著令 $m := m_x + P_{xy} P_y^{-1}(y-m_y)$ 可得
\[{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\]注意到
\[\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right] = \det {P_y}\det P\]故可得
\[{p_{x|y}}\left( {x|y} \right) = \frac{1}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det P} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right) = n(x,m,P) \ \ \ \ \ \square
\]
如果要推導 最佳估測器 我們需要上述結果衍生：
================
FACT 1': Joint independent normals
若 $p_{x|z} (x|z) = n(x, m_x, P_x)$ 為 normal ，令 $y \sim N(m_y, P_y)$ 且 與 $x,z$ 彼此獨立 則 conditional joint density of $(x,y)$ given $z$ 為
\[{p_{x,y|z}}\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]|z} \right) = n\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&0\\
0&{{P_y}}
\end{array}} \right]} \right)\]================
================
FACT 2': Linear Transformation of a Normal
若 $p_{x|z}(x|z)= n(x,m, P)$ 且 $y$ 為 $x$ 的線性轉換；亦即 $y = Ax$ 則
\[{p_{y|z}}(y|z) = n(y,Am,AP{A^T})\]================
================
Preliminary
回憶若 $x$ 為 mean $m$ 且 variance $\sigma^2$ 的 Normal 隨機變數 則 我們表示為 $x \sim N(m, \sigma^2)$ 。其 probability density function 可記做
\[
\frac{1}{\sqrt{2 \pi} \sigma} \exp(\frac{-1}{2 \sigma^2}(x-m)^2)
\]現在若我們拓展上述結果到 多個隨機變數 (又稱 random vector 隨機向量) 的情況，令 $x$ 為Normal 隨機向量 記做
\[
x \sim N(m,P); \;\;
p_x(x) := n(x,m,P)
\]上述符號 表示 $x$ 為 normal distributed 且 mean vector $m$ 與 convariance matrix $P$。另外 $n(x,m,P)$ 表示 normal probability density function
\[\;n\left( {x,m,P} \right): = \frac{1}{{{{\left( {2\pi } \right)}^{n/2}}{{\left( {\det P} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right)\]
Comment:
若 $x \in \mathbb{R}^n$ 則 mean vector $m \in \mathbb{R}^n$ 且 convariance matrix $P \in \mathbb{R}^{n \times n}$ 且為 實數 對稱 正定 矩陣 (正定條件用以確保 $P^{-1}$ 存在，使得上述的 probability density function 可以被定義)。若 $P$ 不為對 正定 我們稱為 singular normal distribution 或稱 degenerate normal。
Example
若 $n=2$ 則我們可以繪製 normal density function； 比如說
\[m = \left[ {\begin{array}{*{20}{c}}
0\\
0
\end{array}} \right];\begin{array}{*{20}{c}}
{}&{}
\end{array}{P^{ - 1}} = \left[ {\begin{array}{*{20}{c}}
{3.5}&{2.5}\\
{2.5}&{4.0}
\end{array}} \right]\]則我們可繪製
現在我們看幾個 之後會使用到的基本結果：
================
FACT 1: Joint independent normals
若隨機向量 $x \sim N(m_x,P_x)$ 與 $y \sim N(m_y,P_y)$ 為 normally distributed 且 彼此互為獨立，則其 joint density $p_{x,y}(x,y)$ 如下
\[
p_{x,y}(x,y) = p_{x}(x)p_{y}(y)=n(x,m_x,P_x) \cdot n(y,m_y,P_y)
\]且
\[\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right] \sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&0\\
0&{{P_y}}
\end{array}} \right]} \right)\]================
================
FACT 2: Linear Transformation of a normal
若 $x \sim N(m, P)$ 且 $y$ 為 $x$ 的線性轉換；亦即對任意矩陣 $A$， $y = Ax$ 則
\[
y \sim N(Am , APA^T)
\]================
================
FACT 3: Conditional of a joint normal
若 $x,y$ 為 jointly normal distributed (no independent assumption)
\[\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]\sim N\left( {\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)\]則 conditional density of $x$ given $y$ 仍為 Normal 亦即
\[
(x|y) \sim N(m,P)
\]其 probability density function 為 $ p_{x|y} (x|y) = n(x,m,P)$ 其中 conditional mean vector $m$ 與 conditional convariance matrix $P$ 分別為
\[\begin{array}{l}
m = {m_x} + {P_{xy}}P_y^{ - 1}(y - {m_y})\\
P = {P_x} - {P_{xy}}P_y^{ - 1}{P_{yx}}
\end{array}\]================
注意到上式中 conditional mean $m$ 為 random vector (depends on $y$) 。
Proof:
回憶 conditional density of $x$ given y 定義
\[
p_{x|y}(x,y) := \frac{p_{x,y}(x,y)}{p_y(y)}
\]注意到 $(x,y)$ 為 joint normal 故我們有
\[\small \begin{array}{l}
{p_y}\left( y \right): = \frac{1}{{{{\left( {2\pi } \right)}^{{n_y}/2}}{{\left( {\det {P_y}} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)\\
{p_{x,y}}\left( {x,y} \right): = \frac{1}{{{{\left( {2\pi } \right)}^{\left( {{n_x} + {n_y}} \right)/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]} \right)
\end{array}\]亦即
\[\begin{array}{l}
{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \frac{{\frac{1}{{{{\left( {2\pi } \right)}^{\left( {{n_x} + {n_y}} \right)/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]} \right)}}{{\frac{1}{{{{\left( {2\pi } \right)}^{{n_y}/2}}{{\left( {\det {P_y}} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right]}^T}{{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]}^{ - 1}}\left[ {\begin{array}{*{20}{c}}
{x - {m_x}}\\
{y - {m_y}}
\end{array}} \right] - {{\left( {y - {m_y}} \right)}^T}{P_y}^{ - 1}\left( {y - {m_y}} \right)} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}
\end{array}\]注意到若我們取 $P:= P_x - P_{xy}P_y^{-1}P_{yx}$ 則 利用 Matrix inversion Lemma 可知
\[{\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]^{ - 1}} = \left[ {\begin{array}{*{20}{c}}
{{P^{ - 1}}}&{ - {P^{ - 1}}{P_{xy}}P_y^{ - 1}}\\
{ - P_y^{ - 1}{P_{yx}}{P^{ - 1}}}&{P_y^{ - 1} + P_y^{ - 1}{P_{yx}}{P^{ - 1}}{P_{xy}}P_y^{ - 1}}
\end{array}} \right]\]將此結果帶回我們可得
\[\begin{array}{l}
{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( \begin{array}{l}
{\left( {x - {m_x}} \right)^T}{P^{ - 1}}\left( {x - {m_x}} \right) - 2{\left( {y - {m_y}} \right)^T}P_y^{ - 1}{P_{yx}}{P^{ - 1}}\left( {x - {m_x}} \right)\\
\begin{array}{*{20}{c}}
{}&{}
\end{array} + {\left( {y - {m_y}} \right)^T}\left( {P_y^{ - 1}{P_{yx}}{P^{ - 1}}{P_{xy}}P_y^{ - 1}} \right)\left( {y - {m_y}} \right)
\end{array} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {{{\left( {x - {m_x}} \right)}^T} - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]}^T}{P^{ - 1}}\left[ {{{\left( {x - {m_x}} \right)}^T} - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\\
\begin{array}{*{20}{c}}
{}&{}&{}
\end{array} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}\left( {{{\left[ {\left( {x - {m_x}} \right) - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]}^T}{P^{ - 1}}\left[ {\left( {x - {m_x}} \right) - {P_{xy}}P_y^{ - 1}\left( {y - {m_y}} \right)} \right]} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}
\end{array}\]接著令 $m := m_x + P_{xy} P_y^{-1}(y-m_y)$ 可得
\[{p_{x|y}}\left( {x|y} \right) = \frac{{{p_{x,y}}\left( {x,y} \right)}}{{{p_y}\left( y \right)}} = \frac{{{{\left( {\det {P_y}} \right)}^{1/2}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right)}}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)}^{1/2}}}}\]注意到
\[\det \left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right] = \det {P_y}\det P\]故可得
\[{p_{x|y}}\left( {x|y} \right) = \frac{1}{{{{\left( {2\pi } \right)}^{{n_x}/2}}{{\left( {\det P} \right)}^{1/2}}}}\exp \left( { - \frac{1}{2}{{\left( {x - m} \right)}^T}{P^{ - 1}}\left( {x - m} \right)} \right) = n(x,m,P) \ \ \ \ \ \square
\]
如果要推導 最佳估測器 我們需要上述結果衍生：
================
FACT 1': Joint independent normals
若 $p_{x|z} (x|z) = n(x, m_x, P_x)$ 為 normal ，令 $y \sim N(m_y, P_y)$ 且 與 $x,z$ 彼此獨立 則 conditional joint density of $(x,y)$ given $z$ 為
\[{p_{x,y|z}}\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]|z} \right) = n\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&0\\
0&{{P_y}}
\end{array}} \right]} \right)\]================
================
FACT 2': Linear Transformation of a Normal
若 $p_{x|z}(x|z)= n(x,m, P)$ 且 $y$ 為 $x$ 的線性轉換；亦即 $y = Ax$ 則
\[{p_{y|z}}(y|z) = n(y,Am,AP{A^T})\]================
================
FACT 3': Conditional of a joint normal
若 $x,y$ 為 jointly normal distributed (no independent assumption)
\[{p_{x,y|z}}\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]|z} \right) = n\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)\]則 conditional density of $x$ given $y,z$ 仍為 Normal ，記為
\[
p_{x|y,z} (x|y,z) = n(x,m,P)
\]其中 conditional mean vector $m$ 與 conditional convariance matrix $P$ 分別為
\[\begin{array}{l}
m = {m_x} + {P_{xy}}P_y^{ - 1}(y - {m_y})\\
P = {P_x} - {P_{xy}}P_y^{ - 1}{P_{yx}}
\end{array}\]
================
Proof:
由 conditional density of $x$ given $y,z$ 定義可知
\[{p_{x|y,z}}\left( {x|y,z} \right) = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\]現在對等號右方 分子分母同乘 $p(z)$ 可得
\[\begin{array}{l} {p_{x|y,z}}\left( {x|y,z} \right) = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\\ \begin{array}{*{20}{c}} {}&{}&{} \end{array} = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_z}\left( z \right)}}\frac{{{p_z}\left( z \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\\ \begin{array}{*{20}{c}} {}&{}&{} \end{array} = {p_{x,y|z}}\left( {x,y|z} \right) \cdot \frac{1}{{{p_{y|z}}\left( {y|z} \right)}} \end{array}\]則由先前 FACT 3 可計算 $p(y,z)$ 並且帶入 $p(x,y|z)$ 即可求得所求。$\square$
ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design".
若 $x,y$ 為 jointly normal distributed (no independent assumption)
\[{p_{x,y|z}}\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right]|z} \right) = n\left( {\left[ {\begin{array}{*{20}{c}}
x\\
y
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{m_x}}\\
{{m_y}}
\end{array}} \right],\left[ {\begin{array}{*{20}{c}}
{{P_x}}&{{P_{xy}}}\\
{{P_{yx}}}&{{P_y}}
\end{array}} \right]} \right)\]則 conditional density of $x$ given $y,z$ 仍為 Normal ，記為
\[
p_{x|y,z} (x|y,z) = n(x,m,P)
\]其中 conditional mean vector $m$ 與 conditional convariance matrix $P$ 分別為
\[\begin{array}{l}
m = {m_x} + {P_{xy}}P_y^{ - 1}(y - {m_y})\\
P = {P_x} - {P_{xy}}P_y^{ - 1}{P_{yx}}
\end{array}\]
================
Proof:
由 conditional density of $x$ given $y,z$ 定義可知
\[{p_{x|y,z}}\left( {x|y,z} \right) = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\]現在對等號右方 分子分母同乘 $p(z)$ 可得
\[\begin{array}{l} {p_{x|y,z}}\left( {x|y,z} \right) = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\\ \begin{array}{*{20}{c}} {}&{}&{} \end{array} = \frac{{{p_{x,y,z}}\left( {x,y,z} \right)}}{{{p_z}\left( z \right)}}\frac{{{p_z}\left( z \right)}}{{{p_{y,z}}\left( {y,z} \right)}}\\ \begin{array}{*{20}{c}} {}&{}&{} \end{array} = {p_{x,y|z}}\left( {x,y|z} \right) \cdot \frac{1}{{{p_{y|z}}\left( {y|z} \right)}} \end{array}\]則由先前 FACT 3 可計算 $p(y,z)$ 並且帶入 $p(x,y|z)$ 即可求得所求。$\square$
ref: J. B. Rawlings and D. Q. Mayne, "Model Predictive Control: Theory and Design".
訂閱：
文章 (Atom)
[Claude] 國小數學加減乘除法計算小遊戲：數學怪獸大亂鬥
心血來潮用 Anthropic Claude Opus 4.6 做的簡單國小數學乘除法計算小遊戲，感嘆AI工具之強大與便利。原本可能要耗時幾天的工作轉眼就完成，時代的巨輪確實在飛速轉動。 數學怪獸大亂鬥(Math Monster Brawl)對戰的國小數學 加減乘除 小遊戲連結...
-
這次要介紹的是數學上一個重要的概念： Norm： 一般翻譯成 範數 (在英語中 norm 有規範的意思，比如我們說normalization就是把某種東西/物品/事件 做 正規化，也就是加上規範使其正常化)，不過個人認為其實翻譯成 範數 也是看不懂的...這邊建議把 ...
-
數學上的 if and only if ( 此文不討論邏輯學中的 if and only if，只討論數學上的 if and only if。) 中文翻譯叫做 若且唯若 (or 當且僅當) ， 記得當初剛接觸這個詞彙的時候，我是完全不明白到底是甚麼意思，查了翻譯也...
-
半導體中的電流是由電子(electron)及電洞(hole)兩種載子(carrier)移動所產生 載子移動的方式： 擴散(diffusion) $\Rightarrow$ 擴散電流 (不受外力電場作用) 飄移(drift) $\Rightarrow$ 飄移電流 (...