初級Mathマニアの寝言

数学は色々なところで応用可能であり、多くの人が数学の抽象的な概念の意味や意義を鮮明に知ることができれば今まで以上に面白い物や仕組みが生まれるかもしれません。このブログは数学を専門にしない人のために抽象的な概念の意味や意義を分かりやすく説明することを目的としています。数学を使って何かしたい人のお役に立てたら幸いです。

可制御可観測な線形システム全体の集合は多様体になる

この記事では、可制御・可観測な線形システム全体の集合は多様体になることを説明します。ここで言う線形システムとは、

で解説した

f:id:ogyahogya:20180511094841p:plain

のことですが、状態方程式表現は行列の三つ組  (A,B,C)\in {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} で決定するので、一つの線形システムは  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の一点である\begin{align} (A,B,C)\end{align}

のことだ考えられます。そうすると、線形システム全体の集合とは\begin{align} {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n}\end{align}

のこととなります。この線形システム全体の集合の中には、

で解説した可制御な線形システムや可観測な線形システムが含まれることになり、上の解説記事を注意深く読むと可制御な線形システム全体の集合とは、 \begin{align} S_c:=\{ (A,B,C)\, |\, {\rm rank} \begin{pmatrix} B & AB & \cdots & A^{n-1}B \end{pmatrix}=n\} \end{align}

のことであり、可観測な線形システム全体の集合とは、 \begin{align} S_o:=\{ (A,B,C)\, |\, {\rm rank} \begin{pmatrix} C \\ CA \\ \vdots \\ CA^{n-1} \end{pmatrix}=n\} \end{align}

のことだということになります。つまり、この記事では、

\begin{align} S_c\cap S_o\end{align}

という可制御可観測な線形システム全体の集合が多様体になるということを説明します。多様体については、

を参照してください。

多様体の開集合は多様体になる。

結論から先に言うと、可制御可観測な線形システム全体の集合  S_c\cap S_o が多様体であることは、  S_c\cap S_o が多様体である  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の開集合であることから分かります。そこで、多様体の開集合は多様体になるという事実を解説しておきます。あとで、  S_c\cap S_o {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の開集合だということを説明します。

そもそも開集合とは何かというと、

で説明したように、位相空間の要素のことでした。位相空間  (X,\mathcal{O}) が与えられると、 X の任意の部分集合  Y に対して、

\begin{align} \mathcal{O}_Y:= \{ U\cap Y\,|\, U\in \mathcal{O} \} \end{align}

とすることで、 \mathcal{O}_Y Y の位相になります。このような  \mathcal{O}_Y X の位相  \mathcal{O} から導かれた  Y相対位相と呼んで、位相空間  (Y,\mathcal{O}_Y) (X,\mathcal{O})部分空間と呼びます。多様体は上のリーマン多様体の記事で定義したように任意の異なる2点を開集合で分離できるという性質を持ったハウスドルフ空間でした。実は、ハウスドルフ空間  X の部分空間  Y もハウスドルフ空間になります。なぜなら、任意の  a, b\in Y,  a\neq b に対して、 X はハウスドルフ空間なので、 a\in U, b\in V, U\cap V=\emptyset となる  X の開集合  U, V が存在して、 U':= U\cap Y, V':=V\cap Y Y の開集合で  a\in U', b\in V', U'\cap V'=\emptyset となるからです。

上のことから、多様体  (X,\mathcal{O}) の開集合  Y \mathcal{O} から導かれた相対位相を導入すると開集合  Y はハウスドルフ空間となります。さらに、多様体  X C^{\infty} 級の座標近傍系が  \{ (U_{\lambda}, \phi_{\lambda})\}_{\lambda\in \Lambda} だとすると、これを  Y に制限した  \{ (U_{\lambda}\cap Y, \phi_{\lambda}|U_{\lambda}\cap Y)\}_{\lambda\in \Lambda} Y C^{\infty} 級の座標近傍系になります。したがって、多様体  X の任意の開集合は  X の位相から導かれた相対位相を導入することで多様体になります

距離空間は位相空間になる

応用上もっとも重要な空間は距離が定められた空間であり、その距離を使って位相を導入できる、つまり、空間の異なる点たちの近さを議論できるようになります。あとで線形システム全体の集合にも距離を導入して位相空間としますので、ここでは簡単に距離空間について説明しときます。

距離空間の定義は以下の通りです。

f:id:ogyahogya:20180513132711p:plain

ベクトル空間にノルムが導入されたノルム空間は距離空間になります。ここで、ノルム空間の定義は以下の通りです。

f:id:ogyahogya:20180513132822p:plain

実際に、ベクトル空間  X にノルム  ||\cdot || が導入されているとすると、 X 上の任意の2点  x,y に対して

\begin{align} d_{\rm norm}(x,y):= ||x-y|| \end{align}

と定めると、 X 上に距離  d_{\rm norm}:X\times X\rightarrow {\bf R} が導入できます。この  d_{\rm norm} をノルム  ||\cdot|| から定められる距離と言います

距離はノルムを使わなくても定義できます。例えば、集合  X の任意の2点  x,y に対して

\begin{align} d(x,y) := \begin{cases} 1\quad (x\neq y\, のとき) \\ 0 \quad (x=y\, のとき) \end{cases} \end{align}

と定めると、 X 上に距離  d が導入できますが、これはあるノルムから定められた距離ではありません。 というのは、ノルムから定められる距離  d_{\rm norm} はノルムの定義の(3)から  d_{\rm norm}(\lambda x, \lambda y) = |\lambda| d_{\rm norm}(x,y) を満たしますが、上の距離  d はこれを満たさないからです。

距離空間  (X,d) に対しては次のように距離  d を用いて位相を導入できます。

f:id:ogyahogya:20180513140632p:plain

f:id:ogyahogya:20180513140643p:plain

f:id:ogyahogya:20180513140820p:plain

工学系の研究論文は、距離空間を扱うことが多く、位相について何も言わなくても暗黙のうちに上の意味での距離  d から定められる位相  \mathcal{O}_d が導入されていることがほとんどです。この記事でも距離空間の位相は上のように距離  d から定義された位相を導入するとします。

線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} は多様体になる

線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} が多様体になることを確認しておきます。

まず、線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n}  は次の和(+)とスカラー倍( \cdot)でベクトル空間になります。

\begin{align}  (A_1,B_1,C_1)+(A_2,B_2,C_2) &:=(A_1+A_2,B_1+B_2,C_1+C_2)  \\  \lambda \cdot (A,B,C)&:= (\lambda A,\lambda B, \lambda C) \end{align}

次に、ベクトル空間である線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} に以下のノルムを導入します。任意の  (A,B,C)\in {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} に対して、

\begin{align} ||(A,B,C)||:= ||A||_F +||B||_F + ||C||_F \end{align}

この  ||\cdot ||_Fフロベニウスノルムと呼ばれており、任意のサイズの行列  M に対して

\begin{align} ||M||_F:= \sqrt{\sum_{i,j} m_{ij}^2} \end{align}

というように各成分の2乗和の平方根で定義され、行列の集合のノルムとしてよく使われています。上で説明したように、線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} には

\begin{align} d( (A_1,B_1,C_1), (A_2,B_2,C_2) ) &:= ||(A_1,B_1,C_1)- (A_2,B_2,C_2)|| \\ &=||A_1-A_2||_F + ||B_1-B_2||_F+||C_1-C_2||_F\end{align}

というように上のノルムを用いて距離が定義できて距離空間となります。したがって、上で説明したように、この距離  d を用いて線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} には位相  \mathcal{O}_d が定義でき位相空間になります。

この位相空間  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} がハウスドルフ空間 (任意の異なる2点が開集合で分離できる位相空間) になることは任意の  p_1:=(A_1,B_1,C_1), p_2:=(A_2,B_2,C_2)\in {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} とすると、 B(p_1; d(p_1,p_2)/2), B(p_2; d(p_1,p_2)/2) \in \mathcal{O}_d

\begin{align} B(p_1; d(p_1,p_2)/2)\cap B(p_2; d(p_1,p_2)/2) =\emptyset \end{align}

となることから分かります。さらに、 {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の座標近傍は {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} 自身と行列の3つ組をベクトル化することと同値な微分同相写像  \psi: {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n}\rightarrow {\bf R}^{n(n+m+p)} の組だけの集合  \{({\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n}, \psi)\} であり、これは  C^{\infty} 級の座標近傍系なので線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} は多様体になります。

可制御可観測な線形システム全体の集合  S_c\cap S_o は多様体になる

可制御可観測な線形システム全体の集合  S_c\cap S_o が多様体になることは  S_c S_o が線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の開集合であることから従います。なぜなら、 S_c S_o が多様体  {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n} の開集合なら位相空間の定義から  S_c\cap S_o も開集合となり、多様体の開集合は上で説明したように多様体になるからです。

可制御な線形システム全体の集合  S_c が開集合であることを示します。可観測な線形システム全体の集合  S_o が開集合であることも同様に示せます。 S_c が開集合であるとは、任意の  (A,B,C)\in S_c とすると、ある \delta > 0 が存在して、

\begin{align} B( (A,B,C); \delta ) \subset S_c \end{align}

となることを意味しているので、これを示します。 (A,B,C)\in S_c だとすると

\begin{align} \begin{pmatrix} B & AB & \cdots & A^{n-1} B \end{pmatrix} \end{align}

の行ベクトル  v_1, v_2, \ldots, v_n\in {\bf R}^{1\times nm} は1次独立です。つまり、

\begin{align} (\lambda_1, \lambda_2, \ldots, \lambda_n)\neq (0, 0, \ldots, 0)\quad \Rightarrow \quad  \lambda_1 v_1+\lambda_2v_2+\cdots +\lambda_n v_n \neq 0 \end{align}

が成り立ちます。行ベクトル  v_1,v_2\ldots, v_n たちは行列  A, B の各成分の和と積から構成されているので、

\begin{align}  f_i: {\bf R}^{n\times n}\times {\bf R}^{n\times m}\times {\bf R}^{p\times n}\rightarrow {\bf R}^{1\times nm}\quad {\rm s.t.}\quad v_i= f_i(A,B,C) \end{align}

となる連続関数  f_1,f_2,\ldots, f_n が存在します。連続関数のスカラー倍や連続関数同士の和も連続関数となるので、

\begin{align} f(A,B,C):= \lambda_1 f_1(A,B,C) + \lambda_2 f_2(A,B,C)+\cdots +\lambda_n f_n(A,B,C) \end{align}

で定義された関数  f も連続関数となります。よって、任意の  \epsilon>0 に対して、ある  \delta>0 が存在して、 \begin{align} d( (A,B,C), (A',B',C') ) < \delta \Rightarrow \| f(A,B,C)-f(A',B',C') \|_F < \epsilon \end{align}

となります。これより、 \epsilon を適切に選ぶことで

\begin{align} \begin{pmatrix} B' & A'B' & \cdots & A'^{n-1} B' \end{pmatrix} \end{align}

の行ベクトルが1次独立にできることを示せれば、ある \delta > 0 が存在して、

\begin{align} B( (A,B,C); \delta ) \subset S_c \end{align}

となることが示せたことになります。いま  v_1,v_2,\ldots, v_n が1次独立だったことに注意すると、

\begin{align} \exists \alpha>0\quad {\rm s.t.}\quad ||f(A,B,C)||_F=\alpha \end{align}

が成り立ちます。さらに、\begin{align} | \|f(A',B',C') \|_F - \|f(A,B,C)\|_F | < \| f(A,B,C)-f(A',B',C') \|_F < \epsilon \end{align}

が成り立つので、\begin{align} \| f(A,B,C)\|_F -\epsilon< \| f(A',B',C') \|_F<\| f(A,B,C)\|_F +\epsilon  \end{align}

が成り立ちます。よって、\begin{align} \alpha -\epsilon< \| f(A',B',C') \|_F<\alpha +\epsilon \end{align} が成り立つので、 \epsilon 0<\epsilon<\alpha が成り立つように選ぶと、\begin{align} \| f(A',B',C')\|_F \neq 0 \end{align} が成り立ちます。これは、関数  f の作り方から、\begin{align} \begin{pmatrix} B' & A'B' & \cdots & A'^{n-1} B' \end{pmatrix} \end{align} の行ベクトルが1次独立であることを意味しているので、証明完了です。

このように可制御可観測な線形システムの全体の集合  S_c\cap S_o は多様体になって、特に開集合であることが分かりました。  S_c\cap S_o は開集合なので、ある一つの可制御可観測な線形システムの十分小さな近傍内のすべての線形システムは可制御可観測だということになります。実は、可制御可観測な線形システムの全体の集合  S_c\cap S_o は線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m} \times {\bf R}^{p\times n} の中で稠密であることも示すことができます。つまり、線形システム全体の集合  {\bf R}^{n\times n}\times {\bf R}^{n\times m} \times {\bf R}^{p\times n} の任意の開集合を  U としたとき、 (S_c\cap S_o)\cap U\neq \emptyset が成り立ちます。このことより、一様にランダムに線形システム  (A,B,C) を生成すると、可制御可観測なシステムが得られる可能性が極めて高いということや、可制御可観測でないシステムに微小な摂動が加わっただけで可制御可観測なシステムになるということが言えます。言い換えると、可制御可観測なシステムであっても、実際上はほとんど可制御可観測とは言えないようなシステムもあるということであり、

で紹介したような指標を使って可制御性や可観測性の大きさを測ることが応用上は重要だということになります。

参考文献

記事を書くにあたって以下の本を参考にしました。

(1) 距離空間とかの部分を参考にしました。

集合・位相入門

集合・位相入門

 

 (2) 第2章の定理11の可制御な線形システム全体の集合は開集合であることの証明を参考にしました。

現代の最適制御理論〈上〉 (1974年) (数学叢書〈23〉)

現代の最適制御理論〈上〉 (1974年) (数学叢書〈23〉)

 

 

 

arXivについて

この記事はarXivについてのメモです。

arXivとは

arXivとは理数系の英語の原稿を無料で読むことができる以下のウェブサイトです。

https://arxiv.org/

arXivは以下のような特徴を持ちます。

1) arXivは査読がないので誰でもすぐに掲載できます。一方で、専門誌は査読があるので投稿から掲載まで1年以上かかったりします。

2) 多くの場合、専門誌に投稿した論文もarXivに投稿できます。

arXivへの投稿は研究成果の優先権の確保になるので、似たようなことをやってる人がいそうな場合はできるだけarXivに投稿した方が良さそうです。最新の成果も速く広まりますし。

arXivに画像入りの論文を投稿するときの注意

arXivには、基本的にtexのソースファイルを投稿することになりますが、画像を入れるときは、画像のファイルも投稿することになります。texや画像のファイルを投稿すると、arXivがpdfを作ってくれますが、texで

\usepackage[dvipdfmx]{graphicx}

使っていたらダメ(arXivが生成したpdfは、本来あるはずの画像が消滅した内容になっている)です。

画像を載せるためには、\usepackage[dvipdfmx]{graphicx}ではなく

\usepackage[pdftex]{graphicx}

を使うと良いです(TeXworksだと、これ使うとエラーになんだけど、arXivでは問題なくpdfを作ってくれる)。また、elsevierが提供しているクラスファイルelsarticle.clsを利用している場合、texの中に\usepackage[dvipdfmx]{graphicx}等の図を載せるための命令を書かなくても図は載るので、tex中に\usepackage[dvipdfmx]{graphicx}等の命令を書かないでarXivにelsarticle.clsもアップロードすれば画像を載せることができます

arXivに投稿する画像の入ったepsファイルのサイズを小さくする方法

arXivにアップロードするファイルは一つあたり400KB、合計で1MBの制限があるらしく、制限に引っかかるためepsファイルをアップロードできないことが多いです。そのため、epsファイルのサイズを小さくする方法を知っていると良いです。小さくする方法はたくさんありますが、一例をあげると、epsファイルを「GSview」で開き、Fileからconvertを選びResolution(解像度)を300に設定すると、見た目はほとんど変わらないで良い感じに小さくなったepsファイルが出来上がります。また、epsファイルをpdfファイルに変換してもサイズが小さくなります。epsファイルのpdfファイルへの変換はオンライン上でもできて、例えば以下でもできます。

convertio.co

arXivに投稿した自分の原稿はこんな感じです。

1) 与えられた行列に最も近いグラフラプラシアンを構成する方法を提案した論文。

[1802.06482] Optimal Graph Laplacian

2) 多数のリーダーがいるマルチエージェントシステムのリーダーの数を減らしたときに、リーダーを入れ替えた方がもとのシステムのパフォーマンスに近くなるか否かを論じた論文。結論は、入れ替えない方が近い。

[1802.06479] Optimal leader selection and demotion in leader-follower multi-agent systems

3) 安定な線形システムの H^2 ノルムの意味では最強のモデル低次元化法を提案した論文。 H^2 ノルムの意味ではこれより良い結果になる方法はない。

[1803.01523] Riemannian optimal model reduction of stable linear systems

4) IEEE Transactions on Automatic Control に掲載された論文。このように掲載済みの論文も出版社がarXivにアップすることを許している場合がある(出版社によっては掲載してから2年たったら許すとかがある)。ここで説明している正定値対称行列全体の集合のリーマン幾何学は制御理論だけでなく、信号処理、情報幾何などの研究でも役立つと思われる。

[1803.04097] Structure-preserving $H^2$ optimal model reduction based on Riemannian trust-region method

5) リーマン多様体上の最適化手法を用いた状態遷移行列に対称性をもつ連続時間線形システムの同定法

[1804.10379] Riemannian optimal system identification method of linear continuous-time systems with symmetry

 

グラフラプラシアン

この記事では、電力網のネットワーク、交通網のネットワーク、人間関係のネットワーク、神経ネットワーク、遺伝子ネットワークのようなネットワークシステムの性質を解析する際に重要なグラフラプラシアンについて解説します(枝に向きのないネットワークだけ解説します)。

グラフ、隣接行列、次数行列

下図のように節点と枝から構成されるネットワークを数学的に表現するには、グラフという概念が役立ちます。

f:id:ogyahogya:20180119164435p:plain

グラフとは、節点の集合  V と枝の集合  E\subset V\times V の組  (V, E) のことです。例えば上のネットワークだと節点集合が \{1, 2, 3, 4, 5, 6\} で枝集合が  \{(1,2), (1,3), (1,5), (2,1),(2,3), (3,1), (3,2), (3,6), (4,5), (4,6), (5,1), (5,4), (6,3),(6,4)\} です。このように  (i,j)\in E なら  (j,i)\in E が成り立つグラフを正確には無向グラフといいます(この記事では、無向グラフだけを説明します)。節点は頂点、枝は辺とも呼ばれます。

上のグラフはそれぞれの枝が同等の重要度を持っているとすると、次のような行列で表現できます。

f:id:ogyahogya:20180123123124p:plain
つまり、上のグラフは枝  (1,2) が存在するので行列の  (1,2) 成分と  (2,1) 成分のところを1として、枝  (1,4) は存在しないので行列の  (1,4) 成分と  (4,1) 成分を0としています。他の成分も同様に定義すると、上のような対称行列が得られます。逆に上の対称行列から上のグラフが構成できます。もしグラフのそれぞれの枝で重要度が異なる(重み付きグラフという)なら、重要度に応じて1以外の正の実数を行列の成分に定義することで重み付きグラフを対称行列によって表現できます。

より数学的には、 グラフ  G=(V,E) の枝と重み(重要度)との対応関係を表す関数  w:E\rightarrow {\bf R}_+ をグラフ G 上の重み関数といいます(  {\bf R}_+ は正の実数の集合)。ただし、 (i,j)\in E なら  w ( (i,j) )=w( (j,i) ) を満たします。この重み付きグラフを  (G,w) と書きます。このような重み付きグラフ  (G,w)を行列によって表現した対称行列  A=(a_{ij})隣接行列といい、次のように定義されます。

f:id:ogyahogya:20180123162009p:plain

通常のグラフ  G はすべての  (i,j)\in E に対して  w( (i,j) )=1 となる重み付きグラフのことです。

重み付きグラフ  (G,w) が与えられたとき、節点  i\in V次数 d_i:= \sum_{j=1}^n a_{ij} で定義します。ただし、 n はそのグラフの節点の数です。次のように各節点の次数を対角上に並べた対角行列を重み付きグラフ  (G,w)次数行列と呼びます。\begin{align}
D={\rm diag}( d_1, d_2,\ldots, d_n)
\end{align}

例えば上のグラフが与えられたとすると、次数行列は  D={\rm diag} (3, 2, 3, 2, 2, 2) です。

グラフラプラシアン、接続行列

重み付きグラフ  (G,w) が与えられたとします。このとき、隣接行列  A と次数行列  D を用いて定義される次の行列  L\in {\bf R}^{n\times n} を重み付きグラフ  (G,w)グラフラプラシアンといいます。\begin{align}
L:=D-A
\end{align} グラフラプラシアンの定義からグラフラプラシアン  L\in {\bf R}^{n\times n} は固有値 0 を少なくとも一つもち、1をすべての成分に並べたベクトル  {\bf 1}\in {\bf R}^n は固有値 0 に対応する固有ベクトルとなります。なぜなら、 L{\bf 1} = 0 となるからです。

グラフラプラシアンは接続行列と呼ばれる行列によっても特徴付けられます。接続行列を定義するために、重み付き無向グラフ  (G,w)枝に一意な番号と任意に向きを付けます。例えば、はじめに与えた例のグラフの枝に次のように番号と向きをつけてみましょう。

f:id:ogyahogya:20180124100158p:plain

 この枝に番号と向きの付いたグラフは次の行列と対応付けることができます。

f:id:ogyahogya:20180124104614p:plain
つまり、枝が節点から出るときに 1 を、枝が節点に入るときに -1 を、どちらでもないときに 0 を割り当てるというルールで枝に番号と向きの付いたグラフを行列表現できます。このような行列を接続行列といい、正確には次のように定義されます。
f:id:ogyahogya:20180124114512p:plain
重み付き無向グラフの接続行列を定義したとき、各枝には番号が付いています。したがって、枝の上に定義されていた重み関数  w:E\rightarrow {\bf R}_+ から各枝に付いている番号の上の関数  w_e w_e = w( (i,j) ) = w( (j,i) ) (  e= (i,j) = (j,i) ) というように定義されます
 
重み付き無向グラフ  (G,w) のグラフラプラシアン  L\in {\bf R}^{n\times n} と接続行列  B\in {\bf R}^{n\times m} は、枝の番号の上に定義された関数  w_e を用いて次の関係で結ばれます。
f:id:ogyahogya:20180124122305p:plain
この等式は枝に任意に向きを付けた接続行列から、グラフラプラシアンという枝の向きに独立な行列を構築できることを意味しています。また、この等式から重み付き無向グラフ  (G,w)グラフラプラシアンは固有値 0 を少なくとも一つ持つ半正定値対称行列であることが即座に分かります。この等式の証明は以下の通りです。

f:id:ogyahogya:20180124164228p:plain

f:id:ogyahogya:20180124164354p:plain

f:id:ogyahogya:20180124164514p:plain

 グラフラプラシアンという名前の由来

グラフラプラシアンという名前はベクトル解析で登場するラプラシアンという作用素に由来します。ここでは重みのない無向グラフ  (G,w) 、つまり任意の  (i,j)\in E に対して  w( (i,j) ) =1 となる無向グラフのグラフラプラシアンがどのようにラプラシアンに関係しているかを説明します。

まず、グラフラプラシアン  L\in {\bf R}^{n\times n} のベクトル  \psi\in {\bf R}^n への作用は次のようになることに注意します。

f:id:ogyahogya:20180125133359p:plain

次に、以下のような向きと重みのない格子グラフを考えます。

この格子グラフの節点を格子座標  (x,y) で表し、節点の集合を  V とすると、上図のようにグラフの境界から離れた節点の次数は4なので、グラフラプラシアン  L の節点集合上の関数  \phi: V\rightarrow {\bf R} への作用は、グラフラプラシアンのベクトルへの作用のアナロジーで、\begin{align}
(L\phi)(x,y) = 4\phi (x,y) -\phi(x+1,y) -\phi(x-1,y)-\phi(x,y+1)-\phi(x,y-1)
\end{align}となると考えられます。
さらに今度は同じ記号を使ってユークリッド空間  {\bf R}^2 上のテイラー展開可能な関数 \phi: {\bf R}^2\rightarrow {\bf R} を考えます。このとき、  h を非常に小さいとして、 \phi(x+h,y) \phi(x,y+h) を点  (x,y)\in {\bf R}^2 の周りでテイラー展開すると、\begin{align}
\phi(x+h,y) &= \phi(x,y) + h\frac{\partial \phi}{\partial x}(x,y) + \frac{h^2}{2!}\frac{\partial^2 \phi}{\partial x^2}(x,y)+\cdots \\ \phi(x,y+h) &= \phi(x,y) + h\frac{\partial \phi}{\partial y}(x,y) + \frac{h^2}{2!}\frac{\partial^2 \phi}{\partial y^2}(x,y)+\cdots
\end{align}となるので、\begin{align}
\frac{\partial^2 \phi}{\partial x^2}(x,y) &\approx \frac{\phi(x+h,y)+\phi(x-h,y)-2\phi(x,y)}{h^2} \\ \frac{\partial^2 \phi}{\partial y^2}(x,y) &\approx \frac{\phi(x,y+h)+\phi(x,y-h)-2\phi(x,y)}{h^2}
\end{align}となります。よって、\begin{align}
\frac{\partial^2 \phi}{\partial x^2}(x,y)+ \frac{\partial^2 \phi}{\partial y^2}(x,y) \approx \frac{\phi(x+h,y)+\phi(x-h,y)+\phi(x,y+h)+\phi(x,y-h)-4\phi(x,y)}{h^2}
\end{align}となりますが、 {\bf R}^2 上のラプラシアン  \Delta \Delta= \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2} なので、\begin{align}
(\Delta \phi)(x,y) \approx \frac{\phi(x+h,y)+\phi(x-h,y)+\phi(x,y+h)+\phi(x,y-h)-4\phi(x,y)}{h^2}
\end{align}となります。したがって、グラフラプラシアン  L はラプラシアン  \Delta の(-1)倍の離散化とみなせます

グラフラプラシアンの名前の由来のもう少し深いところ

上で説明したように、グラフ理論のグラフラプラシアンとベクトル解析のラプラシアンには密接な関係があります。ここでは先ほどと同様に重みのない無向グラフを対象として、もう少し深いつながりを見たいと思います。

ベクトル解析のラプラシアン  \Delta は\begin{align} \Delta = {\rm div} ({\rm grad}) \end{align}のことです。ここで、 {\rm grad} は勾配(gradient)の作用素、 {\rm div} は発散(divergence)の作用素を表します。つまり、ラプラシアンは勾配と発散の合成によって表現されます。

このことと同様なことがグラフラプラシアンでも言えて、接続行列の転置が勾配に対応し、接続行列が発散に対応します。つまり、\begin{align} L = BB^T \end{align} が成り立ちます。この関係式は上で既に証明しています(今は重みなしで考えているので  w_1=w_2=\cdots =w_m=1 であることに注意)。つまり、上で証明したグラフラプラシアンと接続行列の関係式は、ベクトル解析のラプラシアンが勾配と発散によって関係づけられていることから自然な関係式なのだということが分かります。

まず、接続行列の転置  B^T\in {\bf R}^{m\times n} が勾配  {\rm grad} に対応することを示します。  B^T \phi = (\phi_1,\phi_2,\ldots, \phi_n)\in {\bf R}^n に作用させると、第k成分は\begin{align} (B^T\phi)_k = \sum_{i=1}^n B_{ik} \psi_i \end{align} となります。接続行列  B の行は節点、列は枝に対応していたので、 k=(i,j)\in E で節点iから節点jに行く向きが付いてると  B の第k列の第i成分は1、第j成分は-1になります。よって、\begin{align} (B^T\phi)_k = \phi_i-\phi_j \end{align}となり、節点集合上に定義された関数の差、つまり、勾配として  (B^T\phi)_k を見ることができます。さらに、 B^T\phi 自身は枝上に定義された関数として解釈できます。したがって、ベクトル解析の勾配はスカラー場(関数)をベクトル場(関数)にしましたが、グラフ上の勾配(接続行列の転置)は節点集合上の関数を枝集合上の関数にすると解釈できます

次に、接続行列  B\in {\bf R}^{n\times m} が発散  {\rm div} に対応することを示します。 \psi = (\psi_1,\psi_2,\ldots, \psi_m)\in {\bf R}^m として  B を作用させると、\begin{align} (B\psi)_i = \sum_{k=1}^m B_{ik}\psi_k \end{align} となりますが、接続行列の定義から、これは\begin{align} (B\psi)_i = \sum_{節点iから出る枝k} \psi_k -\sum_{節点iに入る枝k}\psi_k \end{align}を意味しています。例えば、以下のような感じです。

f:id:ogyahogya:20180126111510p:plain

これは、ベクトル解析の発散のイメージそのものです。ベクトル解析の発散のイメージについては、

EMANの物理学・電磁気学・ガウスの定理の証明

あたりを参考にしてください。

以上をまとめると、以下の対応関係が得られます。

f:id:ogyahogya:20180125165128p:plain

先ほどまで、勾配、発散、ラプラシアンはユークリッド空間上の関数に作用させていましたが、リーマン多様体上の関数に対しても勾配、発散、ラプラシアンが定義でき、上のような対応関係が得られます。リーマン多様体については、以下を参考にしてください。

ogyahogya.hatenablog.co

参考文献

記事を書くにあたって参考にした文献です。

(1) 重み付きグラフとグラフラプラシアンの定義のところを参考にしました。グラフ理論の初歩的な内容の解説と、初歩的なグラフ理論の制御工学の問題への応用が分かりやすく説明してあります。

マルチエージェントシステムの制御 (システム制御工学シリーズ)

マルチエージェントシステムの制御 (システム制御工学シリーズ)

  • 作者: 東 俊一,永原 正章,石井 秀明,林 直樹,桜間 一徳,畑中 健志
  • 出版社/メーカー: コロナ社
  • 発売日: 2015/08/28
  • メディア: 単行本(ソフトカバー)
  • この商品を含むブログを見る
 

 (2) ベクトル解析の勾配、発散、ラプラシアンがグラフ理論の接続行列の転置、接続行列、グラフラプラシアンにそれぞれ対応しているという解説を参考にしました。この本は電気回路のネットワークの性質をグラフ理論の概念を使って色々解析しています。面白いです。

線形代数とネットワーク

線形代数とネットワーク

 

 (3) リーマン多様体上の関数への勾配、発散、ラプラシアンの作用素が詳しく説明されています。

ラプラシアンの幾何と有限要素法 (朝倉数学大系)

ラプラシアンの幾何と有限要素法 (朝倉数学大系)

 

 

最小平均二乗誤差推定値は条件付き期待値

この記事では、下図のように  y を観測してパラメータ  x を推定しようとしたとき、推定したいパラメータ  x と推定値  \hat{x} の二乗誤差  ||x-\hat{x}||^2 の期待値が最小となる推定値(最小平均二乗誤差推定値 \hat{x}条件付き期待値  E(x|y) で与えられ、この推定値は不偏推定値になることを説明します。f:id:ogyahogya:20171108130548p:plain

ここでは、 x\in {\bf R}^n,  y\in {\bf R}^p,  \hat{x}\in {\bf R}^n として、 ||\cdot|| をユークリッドノルムとします。

推定したいパラメータ  xと観測値  y は確率変数だと考える。そうすると、推定値  \hat{x} も確率変数だと考えるのが自然になる。

この記事の中では、推定したいパラメータ  xと観測値  y は確率変数だと考えます。そのとき、推定値  \hat{x} も確率変数だと考えるのが自然になります。これらの理由を以下で説明します。

推定したいパラメータ  x は全知全能の神様にとっては確定値かもしれませんが、「 x \hat{x} か?」と思っている人にとっては確定値ではなく、確率的にしか決められないものかもしれません。推定したいパラメータ  x は確率的に定まるものだと考えることは、数学的には、 x は確率変数だと考えることに相当します。確率変数は下の記事で紹介しています。

ogyahogya.hatenablog.com

推定したいパラメータ  x が確率変数だと考えると、観測値  y も確率変数だと考えることが自然なことが多いです。実際に、観測値は  y=x+v というように、確率変数である  x に確率変数である雑音  v が加算されたものかもしれず、もしそうなら観測値  y  は確率変数となるためです。他にも  y x の関数として表現される多くの場合に、 y は確率変数ということになります。

推定値  \hat{x} は観測値  y の関数だと考えて、関数  f を用いて  \hat{x} = f(y) という関係があると考えるのが自然です。なぜなら、推定値  \hat{x} は観測値  y を参考にして決定するはずだからです。したがって、観測値  y が確率変数なら  \hat{x} も確率変数だと考えることができます(関数  f が可測関数という連続関数を含む応用上とても広いクラスの関数なら、 y が確率変数のときに  f(y) も確率変数となる)。

なぜ  E(||x-\hat{x}||^2) を最小化するのか

この記事の冒頭で、二乗誤差  ||x-\hat{x}||^2 の「期待値」を最小にする  \hat{x} は条件付き期待値  E(x|y) だということを紹介すると書きましたが、なぜ二乗誤差  ||x-\hat{x}||^2 そのものではなく、その「期待値」を最小にする  \hat{x} を考えるのでしょうか?その理由は、二乗誤差  ||x-\hat{x}||^2 \hat{x}=x のときに最小値0となりますが、ここでは  x を確率変数だと考えるため  \hat{x}=x とすることができない上に、 x が不明なので、どうやって  ||x-\hat{x}||^2 を小さくできるかも分からないためです。そのため、代わりに二乗誤差の期待値  E(||x-\hat{x}||^2) の最小化を考えるのです。

 E(||x-\hat{x}||^2) を最小化するとは、もっと正確にはどういうことなのか

二乗誤差の期待値  E(||x-\hat{x}||^2) は関係式  \hat{x}=f(y) の関数  f が与えられたら、次のように計算できます。

\begin{align*} E(||x-\hat{x}||^2) &= E(||x-f(y)||^2) \\ &= \int_{ {\bf R}^p } \int_{ {\bf R}^n} ||x-f(y)||^2 p(x,y) dx dy \end{align*}

この式を見れば分かるように、二乗誤差の期待値  E(||x-\hat{x}||^2) は同時確率分布  p(x,y) が与えられたら関数  f の関数なので、 E(||x-\hat{x}||^2)

\begin{align*} R(f) := E(||x-\hat{x}||^2) \end{align*}

と書けます。このように、二乗誤差の期待値  E(||x-\hat{x}||^2)を最小化するということは、同時確率分布  p(x,y) が任意に与えられたときに関数  R(f) を最小化する関数  f を求めることを意味しています。

※関数  R は関数  f の関数だということから、汎関数(関数の関数)だということになります。

汎関数  R(f) を最小にする関数  f(y) は条件付き期待値  E(x|y) になる

二乗誤差の期待値  E(||x-\hat{x}||^2) 、つまり汎関数  R(f) の最小化を考えてみましょう。まず、 R(f)

\begin{align*} R(f)  &= \int_{ {\bf R}^p } \int_{ {\bf R}^n} ||x-f(y)||^2 p(x,y) dx dy \\ &=\int_{ {\bf R}^p } R_c(f)  p(y) dy \end{align*}

となることに注意します。ここで、

\begin{align*} R_c(f):= \int_{ {\bf R}^n} ||x-f(y)||^2 p(x|y) dx =E(||x-f(y)||^2 | y) \end{align*}

です。

実は、 R(f) を最小にする  f R_c(f) を最小にする  f は一致します。したがって、 R(f) を最小化する  f の代わりに、 R_c(f) を最小化する  f を求めれば良いということになります。

実際に、 R_c(f) を計算してみると以下のようになります。\begin{align*} R_c(f) &= E(||x-f(y)||^2 | y) \\ &= E(|| (x-E(x|y)) + (E(x|y)-f(y))||^2 | y) \\ &= E(||x-E(x|y)||^2 | y) + ||E(x|y)-f(y)||^2 \\ &\geq E(||x-E(x|y)||^2 | y) \end{align*}

最後の不等式で等号が成り立つのは  f(y) = E(x|y) のときなので、 R_c(f) を最小にする関数  f(y) は条件付き期待値  E(x|y):= \int_{{\bf R}^n} x p(x|y) dx だということが分かります。よって、最小平均二乗誤差推定値は  \hat{x} = E(x|y) となります。

※上の証明の中で突然  E(x|y) が出てきました。これは、最小平均二乗誤差推定値は  \hat{x} = E(x|y) となることを知っていないと思いつかないかもしれません。このことを知らないとして、関数  R_c を最小にする関数  f を求めるためには微分  \frac{\partial R_c}{\partial f}=0 を満たす f を求めれば良さそうです。実際に、このへんの話はカルマンフィルタ関係の本によく載っているのですが、多くの本のなかで微分  \frac{\partial R_c}{\partial f}=0 を満たす f を求めています。が、実はこの微分といっているものは通常の微分だと考えることができません。なぜなら  f は数値ではなく関数だからです。ということで、微分  \frac{\partial R_c}{\partial f} を計算するといっても通常の意味の微分は計算できません。では、どうするかというと汎関数  R_c変分(汎関数微分)を計算し、その変分がゼロになる関数  f が答えだということになります。つまり、力学でオイラー・ラグランジュ方程式を導出する際にラグランジアンの作用積分の変分を計算したように、汎関数  R_c の変分を計算すれば良いということになります。では実際に汎関数  R_c の変分を計算してみましょう。汎関数   R_c の点  f における  h 方向への微分は

\begin{align*} \delta R_c(f)[h] := \lim_{\epsilon\rightarrow 0} \frac{R_c (f+\epsilon h) -R_c(f)}{\epsilon} \end{align*} となります(上記の  \delta R_c(f) がゼロとなる関数  f を見つけることを変分問題  \delta R_c(f)=0 と書く)。 汎関数  R_c(f) の定義から、

\begin{align*} R_c(f+\epsilon h) -R_c(f) = \int_{{\bf R}^n} \left( 2\epsilon (x-f(y))^T h(y) + \epsilon^2 ||h(y)||^2 \right) p(x|y) dx \end{align*}

となります。よって、

\begin{align*} \delta R_c(f)[h] &= 2 \int_{ {\bf R}^n} (x-f(y))^T h(y) p(x|y) dx \\ &= 2(E(x|y)-f(y))^T h(y) \end{align*}

ということになります。したがって、 f(y)=E(x|y) のとき  \delta R_c(f)=0 となることが分かりました。

最小平均二乗誤差推定値  \hat{x} = E(x|y) は不偏推定値

 一般に、推定したいパラメータ  x と推定値  \hat{x} の誤差の期待値  E(x-\hat{x}) がゼロになるとき推定値  \hat{x} は推定したいパラメータ  x の不偏推定値だと言います。つまり、推定したいパラメータの期待値と推定値の期待値が一致しているときに、その推定値は不偏だと言うのです。

最小平均二乗誤差推定値  \hat{x}=E(x|y) は次のように不偏推定値だということが分かります。

\begin{align*} E(\hat{x}) &= E( E(x|y)) \\ &= \int_{{\bf R}^p} \left( \int_{{\bf R}^n} x p(x|y) dx \right) p(y) dy \\ &= \int_{{\bf R}^p} \int_{{\bf R}^n} x p(x,y) dx dy \\ &= \int_{ {\bf R}^n} x \int_{ {\bf R}^p} p(x,y) dy dx\\ &= \int_{{\bf R}^n} x p(x) dx = E(x) \end{align*}

参考文献

記事を書くにあたって参考にした文献です。

 (1) 理論的なことがしっかり書いてる。

応用カルマンフィルタ

応用カルマンフィルタ

 

この記事の中では、二乗誤差  ||x-\hat{x}||^2 の期待値の最小化を考えましたが、この本の中では、他にも絶対誤差の期待値や一様誤差の期待値の最小化も考えられており、条件付き確率分布  p(x|y)ガウス分布に従うなら、絶対誤差と一様誤差の期待値の最小値を与える関数は二乗誤差の期待値の最小値を与える条件付き期待値  E(x|y) に一致することが説明されています。ということで、条件付き期待値  E(x|y) は他の評価指標を使っても最も良い推定値になることがあるということが分かります。

ガウス分布は、以下の記事で説明したように中心極限定理と密接に関連した確率分布で応用上よくでてきます。

ogyahogya.hatenablog.com(2) 上の参考文献と同じ著者の本。汎関数の変分を計算すべきところを普通の微分のように計算しているのは誤りであることに注意。この本のように英語の本でも普通の微分のように計算しているものが多々あります。。。

非線形カルマンフィルタ

非線形カルマンフィルタ

 

 

 

可制御性グラミアンと可観測性グラミアン

対象とする線形システム

f:id:ogyahogya:20170630161131p:plain

の可制御性と可観測性の定義は

ogyahogya.hatenablog.com

で紹介しましたが、この記事では可制御性と可観測性の「大きさ」を定量的に測る手段を紹介します。線形システムについては

ogyahogya.hatenablog.com

で紹介しています。この記事の中では行列 A を安定、つまり、 A のすべての固有値の実部は負であると仮定します。

入力エネルギー最小化問題

まず、線形写像 \Psi_c: L^2(-\infty,0] \rightarrow {\bf R}^n

f:id:ogyahogya:20170630161918p:plain

と定義します。ここで、 L^2(-\infty,0] は区間 (-\infty,0] 上で定義された2乗可積分な関数全体の線形空間です。このような空間については、

ogyahogya.hatenablog.com

でも紹介しています。この \Psi_c(u) は無限の過去 t=-\infty に状態が0であるときの現在 t=0の状態 x(0) を表していると解釈できます。実際に、 \lim_{t\rightarrow -\infty} x(t) =0 とすると、行列 A が安定なので、

f:id:ogyahogya:20170704085346p:plain

となるからです。

 v_1, v_2\in L^2(-\infty,0] とするとき、 \langle v_1, v_2 \rangle_2:= \int_{-\infty}^0 v_1^T(t) v_2(t) dt によって \langle \cdot, \cdot \rangle_2 を定義し、 ||v_1||_2:=\sqrt{\langle v_1, v_1 \rangle_2} によって ||\cdot ||_2 を定義します。

このとき、次の入力エネルギー最小化問題を考えます。

f:id:ogyahogya:20170630164230p:plain

ただし、 ||\cdot|| はユークリッドノルムを表します。この問題は、「状態を原点から単位球面上に到達させることが可能な2乗可積分な入力の中で、エネルギーが最小になるものを求めよ。」ということを意味しています。

可制御性グラミアン

上の入力エネルギー最小化問題と次の可制御性グラミアン

f:id:ogyahogya:20170703161119p:plain

 は密接な関係があります。ここでは、その関係を見るための準備をします。

線形作用素 \Psi_c と行列 W_c の間には次の関係が成り立ちます。

f:id:ogyahogya:20170703161504p:plain

ここで、 \Psi_c^* \Psi_c共役作用素です。共役作用素はこちらで説明しましたので、興味がありましたら参照してください。

ogyahogya.hatenablog.com

上の関係が成り立つことは次のように確認できます。

f:id:ogyahogya:20170703163148p:plain

f:id:ogyahogya:20170703163648p:plain

任意の \xi\in {\bf R}^n に対して、上記の式が成り立つことから W_c=\Psi_c\Psi_c^* が成り立つことが証明できました。

可制御性グラミアンの表式から W_c は半正定値対称行列だということが分かります。実は、線形システムが可制御であることと、 W_c が正定値対称行列であることは等価です。ここでは簡単のために、線形システムは可制御であるとして話を進めます。

直交射影

 入力エネルギー最小化問題の解と可制御性グラミアンの関係を、すっきりと理解するために、直交射影の概念を紹介します。

f:id:ogyahogya:20170704093559p:plain

直交射影には、次の性質があります。

f:id:ogyahogya:20170704100020p:plain

入力エネルギー最小化問題の解と可制御性グラミアンの関係

入力エネルギー最小化問題の解は、次のように可制御性グラミアン W_c を用いて表現できることが分かります。

f:id:ogyahogya:20170703164814p:plain

実は、ある直交射影 P: L^2(-\infty,0] \rightarrow L^2(-\infty,0] が存在し、(※)を満たす任意の u\in L^2(\infty,0] に対して、 Pu=u_{最適} となって Pu も(※)を満たすことが示せます。直交射影の定義の下で示したように、 ||Pu||_2 \leq ||u||_2 なので、上の主張が成り立つということになります。

そのような直交射影 P

f:id:ogyahogya:20170704101524p:plain

で与えられます。 Pu が(※)を満たすことは、

f:id:ogyahogya:20170704101724p:plain

から分かります。さらに、★が直交射影であることも少し計算すれば示せます。

さらに、入力エネルギーと可制御性グラミアンの間には、次の関係が成り立つことが分かります。

f:id:ogyahogya:20170704102951p:plain

なお、 W_c の逆行列 W_c^{-1} は、線形システムを可制御であると考えているので存在します。

 可制御性グラミアンの幾何学的解釈

ここでは、入力 u\in L^2(-\infty,0] ( ||u||_2\leq 1) で到達可能な状態を可制御性グラミアン W_c を用いて特徴付けます。この特徴付けは次の事実に基づきます。

f:id:ogyahogya:20170704111337p:plain

f:id:ogyahogya:20170704111350p:plain

上で示したことから、入力 u\in L^2(-\infty,0] ( ||u||_2\leq 1) で到達可能な状態は W_c^{1/2} x ( ||x||\leq 1 )だということが分かりました。

以上のことから「可制御性の大きさ」、すなわち、「制御のしやすさ」と可制御性グラミアン W_c の関係を調べることができます。このことを見るために、単位エネルギー ||u||_2= 1 の入力 u\in L^2(-\infty,0] で到達可能な状態の集合は

f:id:ogyahogya:20170704122935p:plain

であることに注意しましょう。今、 W_c^{1/2} の固有値を \lambda_1,\lambda_2,\ldots, \lambda_n、対応する正規直交固有ベクトルを v_1,v_2,\ldots,v_n とすると、

f:id:ogyahogya:20170704123212p:plain

となり、 \mathcal{E}_c は固有ベクトル v_1,v_2,\ldots,v_n の方向を主軸として固有値 \lambda_1,\lambda_2,\ldots,\lambda_n を主値とする楕円であることが分かります。イメージ図は下のような感じです。

f:id:ogyahogya:20170704123700p:plain

つまり、単位エネルギー ||u||_2= 1 の入力 u\in L^2(-\infty,0] で原点から最も遠くまで到達可能な方向が可制御性グラミアン W_c の最大固有値に対応する固有ベクトルの方向で、2番目に遠くまで到達可能な方向が2番目に大きい固有値に対応する固有ベクトルの方向、3番目に・・・ということになります。言い換えると、最も制御しやすい方向が可制御性グラミアン W_c の最大固有値に対応する固有ベクトルの方向で、最も制御しにくい方向が W_c の最小固有値に対応する固有ベクトルの方向だということです。

以上から、線系システムが可制御であっても、制御のしやすい方向としにくい方向があり、それらを定量的に調べるためには可制御性グラミアン W_c の固有値と固有ベクトルを調べれば良いということが分かりました。

可観測性グラミアン

次に、観測のしやすさを定量的に評価する可観測性グラミアンについて説明します。

今、現在の状態 x(0) x_0 とします。このとき、 \Psi_o:{\bf R}^n\rightarrow [0,\infty)

f:id:ogyahogya:20170704142459p:plain

で定義します。線形写像 \Psi_o [0,t] 上で入力を u=0 としたときの線形システムの出力値と関係付きます。なぜなら、 y(t) = C\exp (At) x_0 = \Psi_o x_0 が成り立つからです。可制御性グラミアンの時と同様の計算で、 \Psi_o は可観測性グラミアン

f:id:ogyahogya:20170704142729p:plain

f:id:ogyahogya:20170704142927p:plain

の関係があることが分かります。ここでは、 L^2 [0,\infty) 上に内積を \langle f_1, f_2 \rangle_2 := \int_0^{\infty} f_1(t)^Tf_2(t) dt で定め、ノルムを ||f||_2:= \sqrt{\langle f, f \rangle_2} で定めることにします。すると、出力は y=\Psi_o x_0 のため、出力エネルギーが

f:id:ogyahogya:20170704143554p:plain

と可観測性グラミアン W_0 を用いて表せることが分かります。可制御性グラミアンの時と同様に、可観測性グラミアン W_o を用いて定義された集合

f:id:ogyahogya:20170704143756p:plain

は楕円であり、主軸は W_o の固有ベクトルの方向で、主値は固有値の平方根ということになります。これは球面上の状態 x_0 にも「観測のしやすさ 」があり、可観測性グラミアン W_o の最大固有値に対応する固有ベクトルの方向の状態が最も観測しやすく、最小固有値に対応する固有ベクトルの方向の状態が最も観測しにくいという意味になります。(ここで、観測のしやすさを出力エネルギー ||y||_2^2 の大小で解釈しています。)

参考文献

 記事を書くにあたって次の本を参考にしました。

A Course in Robust Control Theory: A Convex Approach (Texts in Applied Mathematics)

A Course in Robust Control Theory: A Convex Approach (Texts in Applied Mathematics)

 

 

行列の指数関数と対数関数

この記事では、行列の指数関数と対数関数について解説します。Lie群とLie環という概念を理解するための準備に相当します。

●行列の内積とノルム

 {\bf R}^{n\times n} をn行n列の実数を成分に持つ行列全体の集合とします。このとき、 {\bf R}^{n\times n} に以下の内積とノルムを導入できます。

f:id:ogyahogya:20170202134934p:plain

任意の A, B\in {\bf R}^{n\times n} に対して  ||AB||\leq ||A||\cdot ||B|| が成り立つことがコーシー・シュワルツの不等式を使うことで確認できます。

●行列の指数関数

実数  x の指数関数  e^x x=0 の周りでテイラー展開すると  e^x = \sum_{k=0}^{\infty} \frac{1}{k!} x^k でした。これと形式的に同じになるように行列の指数関数が次のように定義されます。

f:id:ogyahogya:20170202140703p:plain

右辺の無限級数が収束することは以下のように分かります。

f:id:ogyahogya:20170202140816p:plain

●行列の指数関数の性質

行列の指数関数の性質は以下の結果を利用して明らかにすることができます。

f:id:ogyahogya:20170202144600p:plain

解であることは  \gamma(t)=\exp(tX) を左辺と右辺に代入して等しくなることを見れば良いです。一意性は、微分方程式の解の一意性定理から従います。

指数関数の性質を明らかにする前に、二つの行列が可換であるか否かを判定するカッコ積を導入しておきます。

f:id:ogyahogya:20170202145056p:plain

つまり、 [X,Y]=0 なら  X Y の積は可換です。

次の結果が重要です。

f:id:ogyahogya:20170202151228p:plain

これを応用して、以下の結果が得られます。

f:id:ogyahogya:20170202151802p:plain

これは指数関数の値域は一般線形群  {\rm GL}(n,{\bf R}) であることを意味しています。

f:id:ogyahogya:20170202152621p:plain

さらに、以下の結果も成り立ちます。

f:id:ogyahogya:20170202152800p:plain

これらより、 X\in {\bf R}^{n\times n} に対して、 \gamma(t)=\exp(tX) は加法群  {\bf R} から  {\rm GL}(n,{\bf R}) への微分可能な準同型写像  \gamma:{\bf R}\rightarrow {\rm GL}(n,{\bf R}) を与えていることが分かります。この曲線  \gamma は実は次のように  \dot{\gamma}(0) を用いて表されます。

f:id:ogyahogya:20170202154506p:plain

●行列の対数関数

実数  x対数関数  \log x x=1 の周りでテイラー展開すると  \log x = \sum_{k=1}^{\infty} \frac{(-1)^{k-1}}{k} (x-1)^k でした。これと形式的に同じになるように行列の対数関数が次のように定義されます。

f:id:ogyahogya:20170206140134p:plain

上のように、一般の行列  A\in {\bf R}^{n\times n} に対しては  A単位行列に近いときに対数関数  \log(A) は定義できます。

次の結果は  \exp 0 の近傍と  I_n の近傍の間の微分同相写像であり、逆写像 \log で与えられることを示しています。

f:id:ogyahogya:20170206140629p:plain

この対数関数を利用して、次の公式が得られます。

f:id:ogyahogya:20170206144838p:plain

 上の公式は X Y が可換でないときには指数関数  \exp X \exp Y の積にカッコ積  [ X, Y ] の影響が出てくることを表現しています。

●参考文献

 記事を書くにあたって、以下の本を参考にしました。

 

等質空間

この記事では、等質空間の概念について説明します。等質空間なるものをなぜ紹介するかというと、また別の記事で実数を成分に持つ正定値対称行列全体の集合  {\rm Sym}_+(n, {\bf R})幾何学的に考えたいからです。そのような集合を考えたい理由は、 {\rm Sym}_+(n, {\bf R}) 上での最適化問題が工学の問題を考えていると自然に出てくるからです。実際の問題の例はまた今度書くと思いますが、この記事では等質空間について説明します。

●群

等質空間の定義を理解するためには、群の概念を知っている必要がありますので、群の定義を確認しておきましょう。

f:id:ogyahogya:20170126125302p:plain

実数や整数は馴染みがあると思いますが、実数全体の集合や整数全体の集合は和に関して群になっています(単位元は両方とも0)。しかし、実数全体の集合は積に関して群になっていません。0の逆元が存在しないからです。同様に、整数全体の集合も積に関しては逆元が存在しないので、積に関して群になっていません。

重要な群の例に次の一般線形群があります。

f:id:ogyahogya:20170126133611p:plain

 {\rm GL}(n,{\bf R}) が通常の行列の積に関して群になっていることは次のように確かめられます。結合法則が成り立つことは、一般線形群の要素が行列であり、行列の積は結合法則を満たすことから言えます。次に、単位元の存在は、単位行列単位元になることから言えます。最後に、逆元の存在ですが、任意の A\in {\rm GL}(n,{\bf R}) に対して \det A\neq 0 なので言えます。

  H\subset G部分群であるとは、 H G の演算によって群になることを言います。一般線形群  {\rm GL}(n,{\bf R}) の部分群の重要な例に次の直交群があります。

f:id:ogyahogya:20170126135942p:plain

●群の左剰余類への分解

群の部分群が与えられたら以下のように左剰余類というものを考えることができます。

f:id:ogyahogya:20170130143542p:plain

左剰余類の概念を使って、考えている群の上に以下のように同値関係を導入することができます。

f:id:ogyahogya:20170130143705p:plain

しがたって、同値関係  \sim による商集合を考えることができます。左剰余類を使った同値関係による商集合は次のように表されることが多いです。

f:id:ogyahogya:20170130143931p:plain

商集合という概念はこちらの記事で解説しましたので、興味があったら見てください。

 ogyahogya.hatenablog.com

 

●群の作用

群は集合の要素を変換する役目があります。つまり、 G を群、 X を集合とした時、写像  f: G\times X \rightarrow X が与えられるという形で群は登場することが多いです。この  f(g,x) を単に  g\cdot x と書くことにします ( g\in G, x\in X)。

f:id:ogyahogya:20170126141331p:plain

群の作用でも特に推移的に作用することが大事です。

f:id:ogyahogya:20170126154633p:plain

定義から、群が集合に推移的に作用していれば、集合の一つの元に群のある要素をかけることによって集合の任意の元が得られることを意味しています。

例えば、一般線形群  {\rm GL}(n,{\bf R}) は正定値対称行列全体の集合  {\rm Sym}_+(n,{\bf R}) に次のように推移的に作用します。

f:id:ogyahogya:20170126154646p:plain

したがって、 {\rm Sym}_+(n,{\bf R}) の任意の元は  {\rm Sym}_+(n,{\bf R})単位元である単位行列 {\rm GL}(n,{\bf R}) のある元を作用させること得られます。

●等質空間

この記事で解説したかった等質空間は次のように定義されます。

f:id:ogyahogya:20170130161010p:plain

上で示したことから、正定値対称行列全体の集合  {\rm Sym}_+(n,{\bf R}) は等質空間ということになります。

これから集合が等質空間だったら何が言えるかを見ていきましょう。そのために、まず集合の要素を動かさない群の要素の集まりである等方部分群を定義します。

f:id:ogyahogya:20170130144849p:plain

集合  X に作用する群  G の等方部分群  G_x は名前の通り、 G の部分群になっています。よって、等方部分群  G_x によって群  G 上に同値関係を導入できて、その商集合  G/G_x を考えることができます。実は等質空間  X は商集合  G/G_x と次のように同一視できます。

f:id:ogyahogya:20170130145330p:plain

 この定理は次のように証明できます。

f:id:ogyahogya:20170130152522p:plain

この定理を使うと、以下のように正定値対称行列全体の集合という等質空間を一般線形群の直交群による商集合として考えられることが分かります。

f:id:ogyahogya:20170130161310p:plain

この事実が正定値対称行列全体の集合のリーマン幾何学を考える際に役立ちます(これについてはそのうち書くと思います。)。

●参考文献

 この記事を書く際に以下の本を参考にしました。

代数概論 (数学選書)

代数概論 (数学選書)