ベクトル値関数の最適化

ベクトル値を出力に持つ多変数関数の最適化をまとめました!あってるかな...

ベクトル値関数を、

f(x) --- (1)

とします。この長さをゼロにすることが、最適化の目標です。

式(1)をテイラー展開します。

f(xx) = f(x) + JΔx --- (2)

ここで、Jはヤコビアンで、以下です。

J = (∂f/∂x) --- (3)

エラーEを、以下で定義します。これがゼロに近づけばよろしい。

E = f(xx)Tf(xx) --- (4)

式(4)に式(2)を代入すると、

E = (f(x) + JΔx)T(f(x) + JΔx) = ΔxTJTJΔx + 2ΔxTJTf(x) + f(x)Tf(x) --- (5)

式(5)をΔxで微分して、それをゼロとおきます。

JTJΔx + JTf(x) = 0 --- (6)

式(6)はΔxに関する一次方程式系なので、線型代数の方法により、解くことができます。オシマイ!
スポンサーサイト

オイラーの定理

多作のオイラー、彼の定理はたくさんありますが、おそらく一番有名なのが、こちら。

e = conθ+ i sinθ --- (1)

神秘的な式です。θ=πとすると、これも有名な(小説の題材にもなった)、

e + 1 = 0 --- (2)

が導けます。基本的な5つの数が、ひとつの方程式で結ばれているという、驚くべきもの。

さて、式(1)の導出ですが、教科書では、テイラー展開を使うものが多いですね。確かにそのとおりなのですが、あまり直観的な説明ではないですね。数式をいじるとそのようになる、という感じなので、いまいち納得感がありません。もっと直観的にわかる導出って、ないんですかね?そうなると、物理的な導出?

Hidden Figures

映画「ドリーム」を観ました。

これは、アメリカでは2016年に公開されたものです。原題は、"Hidden Figures"。表舞台に出なかった人たちの話、という意味でしょうか。

あらすじや評価などは、ネットにいろいろとありますので省きます。私の興味は、コンピュータ黎明期の、ヒトによる手計算の時代です。要するに、エンジニアリング用途の数値計算。ファインマン自伝の時代と被ります。

最初の場面に出てくる4次方程式の解ですが、2次方程式ふたつに因数分解されていて、解がふたつでてくるので、もうひとつの方程式は、解ナシ(複素解)ということでしょうね、おそらく...(よく見ることができなかった)

最後の場面の、「オイラー法」というのは、これは数値積分のやり方のことだと思います。これも、手計算でやっていたとは...

ジャイロ効果 (2)

広瀬茂男先生著「ロボット工学―機械システムのベクトル解析」、「4.7 ジャイロ効果の再考」、に進みます。改訂版で追加されたところです。以下の複雑な式を使います。

N = I (ω'* + ω0×ω) + ω× --- (1)

角速度ωxで回る円板を、角速度ωyで回すと、z軸周りに回転する、というのがジャイロ効果です。ここで、

ω = ωx + ωy --- (2)
ω0 = ωy --- (3)

とおくのがミソ。これで式(1)を計算すると、z軸周りのトルクが計算できます。

でも、結果だけ得たいのであれば、

N = ωy×x --- (4)

としても、同じ結果が得られます。実際のところ、通常は式(4)で求めると思います。敢えて式(1)を使う理由ですが、ジャイロ効果をより明確に理解するためだそうです(改訂版序文より)。

ジャイロ効果

広瀬茂男先生著「ロボット工学―機械システムのベクトル解析」、改訂版を読んでいます。

http://kanouy.blog9.fc2.com/blog-entry-1864.html

「4.6 伸縮回転する角運動量ベクトルの微分」、のところで、複雑な式が出てきます。すなわち、

N = I (ω'* + ω0×ω) + ω× --- (1)

同書では、式(4.42)です。これは複雑です。通常では、

N = Iω' + ω× --- (2)

程度ですけどね(これも複雑だが)。では、なぜ式(1)なのか?これは、次の「4.7 ジャイロ効果の再考」、で使われます。初版にはなかったところです。

湘南工科大学非常勤講師 (6)

湘南工科大学の非常勤講師、金曜日の2コマめです。この日は、8時にオープンする、湘南T-SITE内、スタバ湘南蔦屋書店、に入り、準備をするのが慣例となりました。楽しい時間です。

前回(2017年10月20日)は雨。以下、スタバ店員さん(かわいい女性)との、注文時の会話。

女性:今日は寒いですね~
私 :雨ですからね...
女性:このあたり(T-SITEのこと)にお住みですか~
私 :いえ、もう少し先です。
女性:どちらに行かれるんですか~
私 :湘南工科大というところです。
女性:えっ、なにを教えられてるんですか~
私 :数学です。
女性:....

やはり、「数学」というのは、雑談ではタブーであることを、再認識いたしました。この単語を出すと、円滑な会話がストップ。変人とも思われますから、気をつけよ!

Levenberg–Marquardt法

Levenberg–Marquardt法(LM法)を使っているのですが、よくわからないところがあります。

まず、ベクトルでのコスト関数を設定します。そして、それをパラメタで偏微分します。偏微分の方向にパラメタを変化させてやればよいはず。

しかし、もしも、パラメタが本来と逆の方向を指していたらどうでしょう?これはあり得ます。なぜかというと、ベクトルでのコスト関数は、符号を反転させてやっても、これの二乗和は変わらないからです。対して、偏微分は符号が反転する。

というわけで、もしも、偏微分の方向でコストが増大すれば、逆方向も探索しないといけないと思うのですが、通常の実装ではどうもそうはしないみたいですね。なぜならば、正しい方向に行っても、コストが増大することがあるからです。その場合は、いまより少しだけ進む。それでも増大すれば、更に少し進む。こうしていくと、最初の状態となり、決して逆方向には進みません。この理屈はわかりました。

う~ん、困ったな...おそらく、私がどこか、勘違いをしているのだと思います。

関数の近似 (2)

関数の近似ですが、やはり参照するのは、Numerical Recipes。私は、第二版はハードカバーにて、第三版はオンラインにて。

第5章が、"Evaluation of Functions"として、このテーマに当てられています。とにかくこれを読めばよろしい。

最初に、Lanczosの、"Applied Analysis"が参考文献として挙げられています。これは1956年に書かれた古典です。Lanczonの本は、"The Variational Principles of Mechanics"を持っていて、これを読めば、力学の全てが解るわけですが、読めていないです。残念...

"Applied Analysis"は購入しようと思います。古いものですが、このあたりの内容は、現代でもさほどの変化はないと思います。

関数の近似

いまさらながら、関数の近似技術が重要であると、再認識中。

たとえば、三角関数、sinθですが、スマホ電卓でも正確に計算できます。でも、

sinθ = x - x3 / 6 --- (1)

と近似してやると、数値的な洞察を得ることができるわけです。エンジニアリングでは、具体的な数値、すなわち数感が大切。

式(1)はテイラー展開ですが、関数近似の一般理論みたいなのはないのでしょうか?数学ではなく、応用技術でしょうかね、このあたりは...

KKT条件

金谷健一先生の「これなら分かる最適化数学」で、KKT条件の復習をしています。

7.2 ラグランジュ乗数、のところです。制約式が等式であれば、通常のラグランジュ未定乗数法を使えばよいのですが、制約式が不等式の場合は、難しくなります。要するに、

λigi(x) = 0 --- (1)

が成り立つということですが、これは解釈が難しいですね。サポートベクタのような、具体的な事例を勉強したほうが、よいかもしれませんね。

湘南工科大学非常勤講師 (5)

湘南工科大学の非常勤講師、2次関数のおさらいをいたしました。

2次関数で重要なのは、平方完成です。これができれば、極値および解の公式が導けます。つまり、2次関数のすべてが理解できる。

数学はここまででよいのでしょうが、実学の大学なので、やはり応用をやらないといけませんね。思案した結果、物理(ニュートン力学)をやろうと思いました。情報系なので、物理の知識はあまり期待できない(?)のですが、ニュートン力学は日常なので、万人の知識として必須と思っています。

平方完成 (4)

2次関数の平方完成は、簡単です。すなわち、

y = ax2 + bx + c --- (1)

において、

y = a(x + (b/2a))2 + c - b2/(4a) --- (2)

では、xがベクトルxとなったとき、式(2)はどうなるのでしょうか?

これを、式(1)のベクトル版からスクラッチでやろうとすると、結構大変です。なので、式(2)を利用します。つまり式(2)は、

y = (x + (1/2)a-1b)a(x + (1/2)a-1b) + c - (1/4)ba-1b --- (3)

と変形できます。こうすれば、あとは、aをA(行列)、bをb(ベクトル)、xをx(ベクトル)に代えます(cはそのまま)。そうすると、

y = (x + (1/2)A-1b)TA(x + (1/2)A-1b) + c - (1/4)bTA-1b --- (4)

となります。転置を入れることに注意。

湘南工科大学非常勤講師 (4)

湘南工科大学の非常勤講師、後学期の初日は、9月22日(金)です。こちらは通常通り、週イチです。自宅から近いのでラク。

後学期の科目は、基礎解析、というものです。初日はガイダンスなので、様子見ですね。このあと、具体的な内容を考えていきたいと思います。シラバスは用意されております。

湘南工科大学非常勤講師 (3)

本日(2017年8月9日)は、湘南工科大学の非常勤講師、夏季集中授業の最終日(五日目)です。科目は線形代数。

台風は通り過ぎ、今日も普通に開講。最後なのでがんばります。90分×4コマ。

湘南工科大学非常勤講師 (2)

本日(2017年8月3日)から一週間(土日除く)、湘南工科大学の非常勤講師をいたします!科目は線形代数です。

この期間も、BLOGは書くと思います。

Lie群

Lie群に再チャレンジ!仕事で必要になってきた。

とある日、数学関連の集まりで、H氏が登場。このおかたは、情報処理分野では重鎮です。なにかの賞も受賞。

話の成り行きで、Lie群の話を私がしたところ、H氏は非常にわかりやすい説明をされました。

後日談としては、H氏は、かなり読者の多いブログを書かれていて、それにLie群の説明を書いてくださいました。これに登場するM氏というのは私のことなので、正しくはK氏です。

exponential mapの謎 (6)

結局のところ、exponential mapから逃れられなくなりました。なぜかというと、ある関数を、exponential mapで微分する必要が生じたからです。

もともとは、PTAM論文で知ったわけですが、ここにもきちんと書かれていないので(5.2節です)、これが参照している下記の本、

V. Varadarajan. Lie Groups, Lie Algebras and Their Representations. Number 102 in Graduate Texts in Mathematics. Springer-Verlag, 1974.

を購入しました。さ~習得するぞ!

平面の当てはめ (2)

湘南工科大学の非常勤講師、夏季集中のカリキュラムを検討中。線形代数です。

内容は、ベクトルと行列です。難しいことはやらず、かつ実用的な題材を提供したい。社会で役立つものです。

そこでですが、まず前半は、三次元の幾何学をやります。ここでの登場人物は、点/平面/直線。ここで、ベクトルを習得してもらう。実際の応用ではn次元になるので、次元数によらない計算技術の習得を目指します。

そして、後半。ここで行列を出すわけですが、ここでの題材をどうするか?いろいろ考えて、共分散行列を登場させようと思いました。難しいって?そうなんですが、実用上、これ以上に重要な行列はないはず。何とか主成分分析まで持っていきたい。なので、固有値問題をかすります。難しいかな?でも、固有値問題って、何に使うのか、知らないひとが多いんです。登場のステージとしては、うってつけだと思います。

ここまで済ませたとして、最後のネタで、前半と後半を融合させたいと考えました。そこで思いついたのが、以下のもの。

http://kanouy.blog9.fc2.com/blog-entry-577.html

これ、締めとして、よくないですか?

重要な行列

線形代数における行列、いろいろと種類がありますが、私見では応用上重要なのは、以下の二つです。

1)対称行列(⊂正方行列)
共分散行列、ヘッセ行列、慣性テンソルなど、エンジニアリングで登場する正方行列のほとんどが対称行列です。対称行列は、直交行列で対角化できるなど、きわめてよい性質がありますので、これをきちんと理解することが重要。

2)非正方行列
これは、最小二乗法で出てくる、いわゆる計画行列(design matrix)。行数はデータ数、列数はパラメタ数。これを基本とした、行列&ベクトル表記による最小二乗法の習得が重要です。

連立一次方程式の解法再構築

線形代数の教科書に記載されている、「連立一次方程式の解法」に異議あり!

通常は、以下のように書かれます。

1)方程式が未知数より多いとき→解なし(不能)
2)方程式と未知数が一致するとき→一意に解が求まる
3)方程式が未知数より少ないとき→解無数(不定)

これはおかしいです。なぜならば、

1)は最小二乗法のことですし、実際にはこの場合が最も多いはずです。実用上極めて重要なものが、解なしとは...。
2)は理論的にはきれいですが、実際にはこのような状況は少ないかも。
3)はいまはやりのスパースです。本によってはℓ2制約を付けるものがあり、これはこれでラグランジュ乗数のよい練習になるのですが、そうして得られた解はスパースではないので、ちょっと不満。

なので、個人的に再構築中です。本も書きますよ!(ウソ)

行列のトレース

行列のトレース(Tr)、使い方がよくわからない。苦手です。定義はもちろん解りますが、応用が利かないのです。

そうしたところに、Bengio+のDL本に、以下の式が書かれてありました。

|A|2 = Tr(AAT) --- (1)

44ページ、同書では式(2.49)です。

式(1)のAは行列ですが、ベクトルにも適用できますから(ベクトルは行列なので)、たとえば、

|x|2 = xTx = Tr(xTx) = Tr(xxT) --- (2)

のように変形ができます。最後の項はシグマを取ると、共分散行列になります。ちょっと見通しがよくなってきた?

テンソルの計算

石井俊全先生著、「一般相対性理論を一歩一歩数式で理解する(2017)」で、テンソルの計算に再挑戦。

キーとなるのは、アインシュタインの縮約(Einstein summation convention)です。慣れると、これは非常に便利です。これを使わなければ、シグマとか、いろいろと面倒くさい記号が入りますが、そういうものは必要なくなります。それから、行列で表現できる二階のテンソルを超えてしまうと、これは行列ではダメなので、テンソルに頼るしかありません。

ただ、縮約で書くと、コンパクトに表現されるものが、実は項数が非常に多かった、ということもありますから、要注意です。

ベクトル場に沿った微分

石井俊全先生著、「一般相対性理論を一歩一歩数式で理解する(2017)」、いかに解りやすいとはいえ、取扱いの内容自体が難しいので、さすがに最終的には、自己責任で理解する必要があります。

たとえば、解りにくいのが、ベクトル場に沿った微分、というヤツ。これは、スカラを微分して作ったヤツ(これは(0, 1)テンソル)と、ベクトル場(これは(1, 0)テンソル)の内積(つまり縮約)なので、(0, 0)テンソル、つまりスカラとなります。

式としてはそうなるのですが、いまいちピンとこないです。方向微分と同じと思えばよいのですが。

ベクトル場の微分、とはまったく違いますから、お間違えなきよう。

クリストッフェルの記号

石井俊全先生著、「一般相対性理論を一歩一歩数式で理解する(2017)」、第5章に、クリストッフェルの記号(Christoffel symbols)の詳細な説明があります。

添え字が上にひとつ、下にふたつあります。合計三つもあるのです。なぜこんなに多いのかというと、この記号は、(直線)座標を曲線座標で2回偏微分したベクトルを、各基底(=1回偏微分したもの)の線形和で表したときの係数だからです。なので、三つなければ係数の識別ができない。

小林昭七先生著、「曲線と曲面の微分幾何(1995)」には、クリストッフェルの記号が出てこないので、ラッキーと思っていたのですが、やはり避けられないようです。幸い、石井先生の本は解りやすいので、今回はチャレンジします。

湘南工科大学非常勤講師

本年度(2017年度)、湘南工科大学の非常勤講師をお引き受けすることになりました。基礎的(かつ実用的)な数学を教えます。

同大学の知人からの依頼です。こちらは、自宅から近いので、なにかとラクであります。

前期と後期があるのですが、前期は、すでにやっている、ものつくり大学と被るので、夏季に集中で行うことにしました。後期は通常の週イチです。

Idempotent Matrix

'Idempotent matrix'という単語を、某書籍で見ました。Idempotent matrixとは、

H2 = H --- (1)

が成り立つ行列Hである、との説明があります。むむ、これは射影(投影)行列のことではないでしょうか。そこで参照したのが、以下の書籍です。

David A. Harville, 'Matrix Algebra From a Statistician's Perspective (2008)'

統計の線型代数本として、評価の高いものです。本書の10章が、Idempotent matrices、というタイトル。ここを読めば、idempotent matrixがわかるのです。まだ読んでいないけれど。6ページのquick courseです。

ラグランジュの未定乗数法 (2)

佐藤一誠先生による、「トピックモデルによる統計的潜在意味解析(2015)」、66ページのあたりで、確率イチの制約条件は実質は必要ない、という驚くべき主張、

http://kanouy.blog9.fc2.com/blog-entry-1703.html

を検証すべく、以前書いた記事、

http://kanouy.blog9.fc2.com/blog-entry-1652.html

に対して、同じ手法を適用してみました。BishopのPRML、第1章のエントロピー、

H = -Σp(xi) log(p(xi)) + λ(Σp(xi) - 1) --- (1)

を最大化する確率密度を求めよ、という問題です。原書では式(1.99)。答えは、

P(xi) = 1/M --- (2)

です(Mはiの取りうる数)。

Hをp(xi)で偏微分して、ゼロとおけばよいのですが(p(xi)をpiと置き換えます)、

∂H/∂pi = -log(pi) - 1 + λ = 0 --- (3)

ここで、λは必要ないということで、式(3)からこれを取り除きます。すると、

∂H/∂pi = -log(pi) - 1 = 0 --- (4)

となります。これをpiで解くと、

pi = 1/e --- (5)

これはどういうことかというと、piはiによらず、定数ということです。式(5)を正規化してやると、式(2)が得られます。というわけで、この問題も、ラグランジュ乗数は必要なかった。ホントですかね...

回転行列のquaternion分解

日本機械学会編「マルチボディダイナミクス(1)-基礎理論」において、オイラー・パラメータというのが登場いたします。第6章。

このオイラー・パラメータというのは、quaternionのことです。分野が違うと、呼び方も違うという、よい見本。ここからは、quaternionと記します。

同書の79ページに、3x3回転行列が、以下のように分解できると記されています。

A = ELT --- (1)

ELは、3x4行列です。これらの成分は、quaternionの成分で構成されます。このように分解すると、いろいろと便利なことがあるというわけです。勉強になります!

ラグランジュの未定乗数法

ラグランジュの未定乗数法(Lagrange multipliers)、理系の教養数学ではあまりきちんと説明されないところですが(大学で習った記憶がない)、応用上、極めて重要であります。

BishopのPRMLには、ラグランジュの未定乗数法の記載がたくさん出てきます(というか、機械学習本では必須の技術。確率の総和がイチであるという制約がありますから)。第1章に、エントロピーの説明があるのですが、そこに、

H = -Σp(xi) log(p(xi)) + λ(Σp(xi) - 1) --- (1)

を最大化するという問題が出てきます。原書では式(1.99)です。PRMLでは、このレベルの解法は、prerequisiteとみなされているので、なんの説明もなく、

P(xi) = 1/M --- (2)

であると記載されます(ただしMはiの取りうる数)。式(2)の意味は、確率が均質なほど、エントロピーは高いということです。これを導出してみましょう。

Hをp(xi)で偏微分して、ゼロとおけばよいので(p(xi)をpiと置き換えます)、

∂H/∂pi = -log(pi) - 1 + λ = 0 --- (3)

となり、piで解くと、

pi = exp(λ-1) --- (4)

となります。常套手段で、両辺でiに関するΣをとると、

1 = M*exp(λ-1) --- (5)

となります。式(4)(5)から、式(2)が求まります。

回転運動の謎 (3)

前回の続きです。

http://kanouy.blog9.fc2.com/blog-entry-1636.html

r' = r'* + ω×r --- (1)

として、

r'' = (r'*)'* + (ω×r)' = r''** + ω×r'* + ωr + ω×r' --- (2)

としましたが、

r'' = (r'*)'* + (ω×r)' = r''** + ω×r'* + (ω×r)'* + ω×ω×r --- (3)

としてもよいはず。こちらのほうが、式(1)を忠実にフォローしていますね。

式(3)の第三項は、

(ω×r)'* = ω'*×r + ω×r'* --- (4)

ですが、

ω'* = ω' --- (5)

ですから、式(3)に式(4)(5)を代入すると、

r'' = r''* + ωr + 2ω×r'* + ω×(ω×r) --- (6)

が得られました。
プロフィール

加納裕(かのうゆたか)

Author:加納裕(かのうゆたか)


[略歴]
1983年3月東京工業大学工学部機械物理工学科卒業
1983年4月(株)図研入社
1987年1月同社退社
1987年2月(株)ソリッドレイ研究所を6名で設立、取締役
1994年3月同社退社
1994年4月(株)スリーディー入社
1996年10月同社取締役
1999年12月上海大学兼務教授
2002年10月同社代表取締役
2009年9月ものつくり大学非常勤講師~現在
2009年10月同社代表退任/退社
2010年1月ソフトキューブ(株)入社~現在(横浜オフィス)
2011年11月甲南大学特別講師
2011年11月関西大学特別講師
2012年11月東京理科大学特別講師
2017年4月湘南工科大学非常勤講師~現在


[業界団体・学会活動]
電気学会・第二期次世代インタラクティブディスプレイ協同研究委員会(幹事)/三次元映像のフォーラム(幹事、監査)/日本バーチャルリアリティ学会・論文委員会(委員)・力触覚の提示と計算研究会(委員)/ACM/SIGGRAPH(Professional Member)/情報処理学会(正会員、CVIM会員)/3Dコンソーシアム(賛助会員)/3DBiz研究会(個人賛助会員)/最先端表現技術利用推進協会(個人会員)/URCF(特別会員)

----------------

前職:/立体映像産業推進協議会(幹事)


[資格]
TOEIC805点
数学検定1級(数理技能)
中型・普自二免許
サッカー4級審判員

最新記事
最新コメント
最新トラックバック
月別アーカイブ
カテゴリ
検索フォーム
RSSリンクの表示
リンク
ブロとも申請フォーム

この人とブロともになる

QRコード
QRコード