最初はちょっと記法覚えるまでのつもりだったけど、どう考えてもあっちのが楽だし、軽いし、というような状況なので、 もう当分はてなダイアリーばかり更新してると思います。

共分散分析

共分散分析を使うと，共変量でもって調整したうえで分散分析をすることができる．たとえば次のようなデータ．

放牧がその後の種子生産に与える影響を調べたデータで，それぞれ初期の根の直径と種子生産量によって散布図が描いてある(統計学:Rを用いた入門書より．データはココ)．

このとき，単純に種子生産量の平均値のみを2つの処理間で比べると，「放牧した方が種子生産量が多くなる」という結果が得られる．

しかし散布図から判断すると，それは放牧区に初期の個体サイズが大きなサンプルが集中していた結果のように見える．たとえば回帰直線を2本引いたら，明らかに放牧区の線は「下」にくるだろう．

それで，2本の回帰直線は違うのかとかそういうことをやるのが共分散分析．

これをやるにはそれぞれの回帰直線が平行という仮定が必要になる．要するにその回帰直線を平均値とし，回帰直線からの偏差を分散として分散分析をするんだから．

なので共分散分析は「平行性の検定(回帰の同質性の検定)」「要因効果の検定(共分散分析)」という2つのステップを踏む．

Rで

Rでは次のように行う．まず平行性の検定．

summary.aov(lm(Fruit~Root*Grazing))

ここでFruitは種子生産量，Rootは根の直径，Grazingはグループを表すベクトル．これにより次のような分散分析表が得られる．

             Df  Sum Sq Mean Sq  F value    Pr(>F)    
Root          1 16795.0 16795.0 359.9681 < 2.2e-16 ***
Grazing       1  5264.4  5264.4 112.8316 1.209e-12 ***
Root:Grazing  1     4.8     4.8   0.1031      0.75    
Residuals    36  1679.6    46.7                 

ここでRoot:Grazingとなっている部分は交互作用を表しており，ここが有意でなければ各回帰直線は平行であるということがいえる．この例ではp値が0.75なので問題ない．

共分散分析は次のように行う．

summary.lm(lm(Fruit~Root + Grazing))

モデル式の書き方がちょっと違う．これにより次のような結果が得られる．

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -127.829      9.664  -13.23 1.35e-15 ***
Root              23.560      1.149   20.51  < 2e-16 ***
GrazingUngrazed   36.103      3.357   10.75 6.11e-13 ***

ここでGrazingUngrazedとなっている部分のEstimateが2つの回帰式がどれだけ離れているかを表しており，ここのp値が有意なら2つの回帰式は有意に離れている(つまり切片が有意に異なる)ということが言える．

(もうちょっと説明をしたいところだけど正直イメージがうまく描けてないのでパスです．青木先生のところなんかをよーく読んだら分かるかもしれません．僕は流し読んだのでわかりません＞＜)

回帰式の多重比較

共分散分析では平行性の検定を行う．で，ここで棄却されなければ次へ進める．

逆に言うと，平行性の検定が棄却されるということは切片はともかくとして2つの回帰式の傾きが有意に異なるということ．

というわけで，群が2つしかない場合にそれらの回帰式の傾きを比べる方法については解決だ．上に書いた方法で平行性を調べたらいい．また，「どこかに差がある」ということが分かればいいんだったらやっぱ同じ方法でいい．

でもやっぱ「どことどこに差があるの？」ってことが気になる．

というわけで，群の中から2つを選び総当たりで平行性の検定を行うことにしよう．そのとき問題になるのは検定の多重性というやつ．

Holmの方法による多重比較

たとえば群が3つあったら対比較は3回やらないといけない．危険率が5％として，すべての群が等しいのにどこかに差がある判定をしてしまう確率は1-0.95^3で0.14くらいになる(実際はちょっと違うらしいけど)．

それで，全体として有意水準5％で検定をしたいときに，対比較の回数に応じて1回1回の対比較における有意水準を低く調整してやることで何とかしようという多重比較法がある．ボンフェローニ(Bonferroni)の方法をはじめとする有意水準調整法というやつだ．

今回はホルム(Holm)の方法というやつにする．これはそれぞれの対比較で得られたp値を小さい順に並べて，i番目のp値に対する有意水準をα/(k-i+1)にするという方法．kは対比較の回数，αは全体としてそれ以下におさめたい有意水準．ボンフェローニの方法だとすべての有意水準をα/kに調整してしまうので過剰に保守的になってしまうけど，ホルムの方法ではもうちょっとゆるい感じになっている．ゆるいといっても全体としての有意水準はα以下に保たれる．

で，やることは決まった．対比較を繰り返して，Holmの方法で有意水準を調整するんだ．

Rで

関数にしてみました．

pairwise.covar <- function(
                           data.f,       #データフレーム
                           g.i=1,        #グループ変数の列
                           i.i=2,        #独立変数の列
                           d.i=3         #従属変数の列
                           )
{
  group <- data.f[[g.i]]                  #グループ変数
  g.name <- levels(group)               #水準名
  independent <- data.f[[i.i]]            #従属変数
  dependent <- data.f[[d.i]]              #独立変数
  g.n <- length(levels(group))          #グループ数
  k <- g.n*(g.n-1)/2                    #組み合わせパターン数
  pvals <- numeric(k)                   #p値用ベクトル
  cmb <- combn(g.name,2)                #全組み合わせ
  cmb.name <- character(k)              #組み合わせの名前（出力用）
  sign <- character(k)                  #アスタリスク用
  for(i in 1:k){
    pvals[i] <- round(anova(
                            lm(dependent[group==cmb[1,i] | group==cmb[2,i]]~
                               independent[group==cmb[1,i]|group==cmb[2,i]]*
                               group[group==cmb[1,i]|group==cmb[2,i]]
                               )
                            )$Pr[3],7)          #p値の計算と格納
    cmb.name[i] <- paste(cmb[1,i],":",cmb[2,i]) #組み合わせ名格納
    if(pvals[i] < 0.001){
      sign[i] <- "***"                  #0.1%以下は***
    }else if(pvals[i] < 0.01){
      sign[i] <- "**"                   #1%以下は**
    }else if(pvals[i] < 0.05){
      sign[i] <- "*"                    #5%以下は*
    }else if(pvals[i] < 0.1){
      sign[i] <- "."                    #10%以下は.
    }else{
      sign[i] <- " "                    #それ以外は" "
    }
  }
  o <- order(pvals)
  cbind(cmb.name[o],
        p.value=(pvals*(k-order(pvals)+1))[o],
        sign[o])                        #結果のベクトル
}

アスタリスクとか飾りだけどさ，ほら，見やすさのために…

使い方は

pairwise.covar(データフレーム, 
				グループ変数の列,
				独立変数の列,
				従属変数の列)

列の指定は省略可．デフォルトでそれぞれ1, 2, 3．

使ってみる．まずデータの作成．

## 仮データの作成
set.seed(5)
x <- 1:20
y1 <- x+rnorm(20)
y2 <- x*2+rnorm(20)
y3 <- x+10+rnorm(20)
group <- as.factor(rep(1:3, rep(20,3)))
independent <- rep(x, 3)
dependent <- c(y1, y2, y3)
mydata <- data.frame(group, independent, dependent)

xが独立変数，y1~3が従属変数．プロットするとこうなる．

凡例忘れた．．．○がy1，△がy2，＋がy3．

見て分かるようにy1とy3は切片が異なるけど傾きは同じ．y2だけ傾きが違う．

これに対してさっきのを使うと

> pairwise.covar(mydata)
             p.value         
[1,] "1 : 2" "0"        "***"
[2,] "2 : 3" "0"        "***"
[3,] "1 : 3" "0.639056" " "  

となる．ちゃんと計算できてる？

ちなみに有意水準を減らす代わりにp値を増やしているので場合によってはp値が1を超えるとかいうことになるかもしれない．あとp値が考えている有意水準より高いからって傾きが等しい＝平行というのにはちょっと無理があるかも．ボンフェローニよりマシとはいえホルムの方法も保守的なのには違いがないので．

Meadowのインストール時にパッケージとしてYaTeXを選択していればこれを.emacsに書き足すだけでいいはず．

もちろん，LaTeX環境は別途準備してある前提で．

拡張子.texのファイルを読み込むとYaTeXが起動し，モードラインに「やてふ」と表示される．

ついでに[C-t l]をdvipdfmxに割り当てた．.dvi作成後に実行で.pdf作成．

Searching for program:no such file or directory, anthy-agent

とか出てauto-autoloads.elなどの読み込みが止まる．止まってもなんか普通にESSとか動くんだけどIMEのON/OFFを切り替えるたびにミニバッファに上の警告が出てエラー音を発する(IMEのON/OFFはちゃんと切り替わる)．

害は今のとこないんだけどなんとかしたかったのでなんとかした．

anthyを切る

IME使うようにしてるのでそもそもいらない．

Meadowのauto-autoloads.elをいくつか止めて快適に - わだいのたけひこのにっき

\Meadow\packages\pkginfo\anthyにあるauto-autoloads.elをリネームして読み込まないようにした．

これでanthy-agantないよ！って怒られなくなった．[-Aあ-\--]とかも出るようになった．起動時に読み込むもの増えた．

でも今度は

error: Invalid value for :char-spec :indian-1-column

とか起動時に言われる．ぐぬぬ…

indian-1-columnをコメントアウト

Meadow起動時のエラー - 教えて!goo

ググってるときに教えて！gooの「一般人」の適当な回答が引っ掛かったときの怒りを表現する言葉が欲しい．

とにかく，indian-1-columnというcharsetがもうないので，読み込まないようにしろと，そういうことらしい．

それでどこに書いてあるんだ，font-setup.elとか検索してもないぞ．と思ったらauto-autload.elの中だった．

Windows Meadow 2.10 多言語環境の設定 - NOX INSOMNIAE

\Meadow\packages\pkginfo\intlfontsの中のauto-autoloads.elを開くと，中に

(indian-1-column ("ind1c16-mule.bdf"));; MuleIndian-1

という記述があるのでコメントアウト．てか消していい気がする．僕はチキンなので消さない．

起動が重くなった

真面目にいろいろ読み込むせいか起動が遅いぞ．いろいろ読み込まないようにしたらいいのだろうか．これから調べる．

*パッケージ減らしたらいくらか早くなった

Emacsのキーバインドを覚える手間を差し引いてもESSを導入する価値はあると思ったのでちょいとメモ．「Rはちょっとくらい使えるけど，Emacsとかよく分からない」というWindowsユーザー向け(というかそもそもMacやLinux使う人はEmacs使えるだろうしESSも自力でインストールできると思うので，RjpwikiのESSの項目なんかを参考にするといいと思います)．

ESSとは？

ESS (Emacs Speaks Statistics) は Emacs や XEmacs から R などの統計解析アプリケーションを便利に使ってしまおう，という Lisp プログラムです．

ESSがどれだけれ能率向上をもたらしてくれるものか，などは「The R Book」p.55--56に叙述されています．(Rjpwiki ESS)

要はEmacsという非常に高機能なテキストエディタがあって，その上でRを使うとき便利なようにいろいろカスタムするプログラムがESSということ．

何を覚える必要があるのか？

1. R
2. Emacs
3. ESS

3つだけ．簡単です．まずは準備を．（ただ最初に書きましたがRの使い方だけはどこかよそで予習しておきましょう．R自体はESSなくてもRguiで動くので．僕に聞くのはいろいろ間違いです＞＜）

R，Meadow，ESSのインストールと設定

.exeクリックで一発！とはいかないのでちょっとがんばる．

Rのインストール

最新版はここにあるはずなので，「Download R 2.8.1 for Windows」などというリンクから.exeを落としてきて実行すればインストールが開始される．

途中の選択肢はESSを使う上ではどう選んでも関係ないはずだけれど，インストールしたフォルダは忘れないように．

Meadowのインストール

Windows用のEmacs実装としてMeadowというものがある．今回はこれにESSを入れて使う．Meadowのページから「ダウンロード」というリンクを辿ると「Netinstall packages」というものが見つかるので，「setup-ja.exe」を落とす．

実行する前に，Cドライブ直下に「Meadow」などの名前でフォルダを作っておき，そこへsetup-ja.exeを移動．

setup-ja.exeを実行してインストール．

• ローカルパッケージディレクトリはC:\Meadow（作成したフォルダ．実行したフォルダへのパスが自動で入るはず）
• インストール先ディレクトリもC:\Meadow（同じく変更の必要はない）
• パッケージはデフォルトでいい気がする

最後にsetup.exeを実行する必要がある．実行すると「どこにある「.emacs」を読むようにする？」などと聞いてくる．画面には[C:\Meadow]などと出ているはず．何も入力せずにEnterを叩けばそこにある.emacsを読むようになる．書き換えてもいい．ただし場所は覚えておく．

これでインストールは終了．デスクトップにショートカットが出現し，それをダブルクリックでMeadowが起動するようになっているはず．

.emacsの作成

Meadowは起動時に.emacsというファイルに書き込まれた設定を読み込む．なのでこれを作らないといけないのだが，「新規作成→テキスト→名前を.emacsに」しても「ファイル名を入力して下さい！」と怒られる．

なのでMeadowから作成する．Emacsの操作が分からない人はまじないの如く以下の作業を行う．

1. Meadowを起動
2. Ctlr+X，Ctlr+Fと連続してタイプ
3. 一番下に「Find file: ~/bin/」と出る
4. 続けて，「Find file: ~/bin/.emacs」と入力してEnter
5. 次に，Ctlr+X，Ctlr+Wと連続してタイプ
6. C:\Meadow\binの中に.emacsというファイルができているはず

完成した.emacsをインストール時に設定した「.emacsを読む場所」に合わせて移動すれば終わり．もちろんわかる人は最初からそこへ作成すればいい．

作ってしまえばこっちのもので，あとはメモ帳で開けば編集できるし，保存もできる．Windowsは難しい．ともかく，Meadowの扱いがわからないうちはメモ帳で編集できるので安心．

ESSのインストール

ESSはたぶんここに最新版があると思う．DownloadからLatest Released Versionって所のzipを落としてきて解凍する．

解凍してできたフォルダをフォルダごと\Meadow\site-lispの中へ放り込む．

.emacsの編集

ESSを使えるようにするため，早速.emacsを編集する必要がある．なんかよくわからなかったらとりあえず次のように書いて保存．

(require 'ess-site)
(setq ess-ask-for-ess-directory nil)
(setq ess-pre-run-hook
	'((lambda ()
	(setq default-process-coding-system '(sjis . sjis))
	)))

(set-language-environment "Japanese")
(setq default-input-method "MW32-IME")
(mw32-ime-initialize)

(defun ess:format-window-1 ()
	(split-window-horizontally)
	(other-window 1)
	(split-window)
	(other-window 1))
(add-hook 'ess-pre-run-hook 'ess:format-window-1)

(setq default-frame-alist
      (append (list '(foreground-color . "azure3")
	'(background-color . "black")
	'(border-color . "black")
	'(mouse-color . "white")
	'(cursor-color . "white")
	)
	default-frame-alist))

（ここのをベースにちょこちょこいじってついでにIMEが動くようにしてあります．）

環境変数の設定

ここまで終わってもまだ動かない．環境変数というものを設定する必要がある．わかっている人はPathへRとMeadowの実行ファイルがあるフォルダ(bin)を，HomeへMeadow\binを追加．以下分からない人へ．

環境変数の設定は，

• XPの場合：コントロールパネル（クラシック表示）→システム→詳細設定タブ など
• Vistaの場合：コントロールパネル（クラシック表示）→システム→左側にある「システムの詳細設定」 など

で見つけられるはず．「ユーザー環境変数」と「システム環境変数」があるけど，ユーザー個人に設定を適用するか，コンピュータ全体に適用するか程度の違いなので以下好きな方へ設定．

1. 変数からPATH(or Path)を選択．なければ「新規」で変数名にPATHと入れる．
2. 変数値にC:\Meadow\bin (RunMW32.exeのあるフォルダ)を追加．新規に作成した場合はそのまま記入．すでに何か書いてある場合，末尾に「;」を記入したのちに続けて記入する．

同様にしてPATHへRの実行ファイルがあるディレクトリ(C:\Program Files\R\R-2.8.1\bin など)を追加し，HOMEへC:\Meadow\bin(.emacsのあるフォルダ)を追加．

確認

拡張子が.RなファイルをMeadowで開くとESSモードで起動．ツールバーに「R」とかのボタンが見えるはず．

M-x R (Alt+x rもしくはEsc x rと順番にキーを叩く)すると，画面が分割されて右下でRが起動．

.Rなファイルが表示されている場所(バッファ)にカーソルがある状態でC-c C-b (Ctrl+c Ctrl+bの順に叩く)すると中身がすべて実行される．結果はRのバッファなどに．

うまくいかなかったら…

どうしよう？よくわかりません＞＜

手順確認→Google先生に聞く　しても分からなかったら聞いてください．でもきっとわかりません．

使い方を覚える

RjpwikiにESSを取り扱ったページがあって， ショートカットなんかは一覧がある(ESSモード)ので参考に練習しましょう．

…

で終わらせたいところだけど，Windowsユーザーにとっての鬼門はMeadow(Emacs)の独特の操作や単語， ショートカットの表記法を覚えるところだと思う．一覧表とか見せられても量がすごいので「覚えられない！無理！」ということになる．

Emacsの操作方法を覚える

Meadow(Emacs)にはチュートリアルが入っている。

ツールバーの[Help]から[Emacs Tutorial]を選ぶとチュートリアルが始まる(ただのテキストだけど)ので， これをやる．

案外長くて丁寧にやると30分くらいかかるけど，覚えるべきことはカーソル操作，消去と挿入，ファイル操作程度なので， 実際それほど大変じゃない．

1回だけじゃ忘れるかもしれないけど，3回もやったら基本的な操作はほとんど覚えられるはず．

「こんな変な操作覚えて大丈夫か？」などという心配があるかもしれないが，大抵のものはMeadowで事足りる， というかMeadowの方が便利なことが多いし，他のエディタがEmacsのキーバインドをサポートしていることも多い(Visual C++とか)．なので心配せず慣れてしまっていいと思う．

ESSの扱いを覚える

操作自体はさっきも書いたRjpwikiのESSのページを参考に覚えたらいい． ショートカットの部分を拡大印刷して手元において置くとかして．

それよりも操作に慣れることが重要．ちょっと慣れると急激にありがたさが分かってくる．

これは実際にいくつかプログラムを書いてみるのが一番いい．

RでProject Eulerなんかいかがですか？

今回、手順としては「PS形式でグラフを作る」→「PSfragでテキストを数式に置換」→「LaTeXで処理」という形をとる。

まー結局はPSfragの使い方なのでRに限った話じゃなくてps/epsファイル全般の話になりますが。

それではお手元にRとLaTeXの環境を用意してご覧下さい。

1．PSfragを準備

まず、PSfragを入手しないといけないのだけど、何か僕の環境にはすでに入っていたのでTeXインストーラーなどを使ってLaTeX一式をまとめてインストールした場合に一緒にインストールされるものに含まれているのだと思う （便利ですねインストーラー！）。もしなかったらCTANなどを探して拾ってきてパスの通った場所に置く。 （この辺正直詳しいことがわからないので聞かないで下さい＞＜）

2．グラフの準備

グラフを準備する。先月作ったデルタ関数のグラフを使いまわそう。

アレは次のようにして描いた。

#余白設定
par(mar=c(3,3,3,1), mgp=c(2,1,0))
#関数定義
Delta <- function(a, x){
 a <- complex(imaginary=1/a)
 Im((pi*(x-a))^-1)
 }
#プロット
curve(Delta(7.5,x), -2, 2, ylim=c(0,10),
main=expression(lim(paste(Im, over(1, pi(x-ia))), a %->% 0)),
xlab="x",
ylab=expression(paste(delta, (x))),
col=1)
par(new=T)
curve(Delta(15,x),  -2, 2, ylim=c(0,10),main="",xlab="",ylab="",col=2)
par(new=T)
curve(Delta(30,x),  -2, 2, ylim=c(0,10),main="",xlab="",ylab="",col=3)
#凡例
legend(
0.5, 10,
box.lty=0,
legend=c(
expression(paste(a=7.5^-1)),
expression(paste(a=15^-1)),
expression(paste(a=30^-1))
),
col=1:3,
lty=c(1,1,1)
)

こんなのが描ける。

…やっぱ数式が汚い。

で、せっかく頑張って作ったこのグラフは使えない。PSfragはPSファイル中のテキストを認識して置換するので、 次のようにしてタイトルや軸ラベルをテキストにしたグラフを作らないといけない。

#プロット
curve(Delta(7.5,x), -2, 2, ylim=c(0,10),
main="main.title",
xlab="xlab",
ylab="ylab",
col=1)
par(new=T)
curve(Delta(15,x),  -2, 2, ylim=c(0,10),main="",xlab="",ylab="",col=2)
par(new=T)
curve(Delta(30,x),  -2, 2, ylim=c(0,10),main="",xlab="",ylab="",col=3)
#凡例
legend(
0.5, 10,
box.lty=0,
legend=c(
"l1",
"l2",
"l3"
),
col=1:3,
lty=c(1,1,1)
)

次のようなグラフが得られる。

ただこれはpng。PSfragが処理できるのはps/epsファイル。 なのでたとえばデバイスをpostscriptにするとか、Windowsのデフォルトのグラフィックデバイスから 「ファイル」→ 「別名で保存」→「postscript」するなどしてepsとして保存しておく。

他に注意として、置換は文字の先頭位置を基準に行われる。なのでたとえばタイトル部分を省略して「t」一文字などにしておくと、 そこを基点に数式が挿入されるのでちょっと右よりになってしまう。

3．LaTeXで処理

パッケージとして呼び出すのはpsfragのほかにgraphicxも必要。psファイルを描く必要があるので。

psfragの使い方は

\psfrag{置換前テキスト}{置換後テキスト}

とこれだけ（オプションも指定できるけど省略）。このコマンドでテキストを置換した後、 psファイルを呼び出してやると置換後の状態で描画される。 置換後のテキスト部分には$$などを用いてLaTeXでいつもやるように数式が書ける。

なので、今の例で行くと次のようなファイルを作る。ここではtexファイルと同じディレクトリに「delta.eps」 という名前で先ほどのpsファイルが保存してある。

\documentclass{jsarticle}
\usepackage[dvips]{graphicx,psfrag}  %パッケージ:graphicx,psfrag
\begin{document}
 \psfrag{main.title}{
$\displaystyle \lim_{a\rightarrow0} \Im \frac{1}{\pi(x-ia)}$
}                                  %テキストの置換
 \psfrag{xlab}{$x$軸}
 \psfrag{ylab}{$\delta (x)$}
 \psfrag{l1}{$7.5^{-1}$}
 \psfrag{l2}{$15^{-1}$}
 \psfrag{l3}{$30^{-1}$}
 \includegraphics[scale=1]{delta}    %epsファイルの呼び出し
\end{document}

で、これをTeXで処理する。dviファイルの段階ではまだ置換は行われず、psファイルまで変換して初めて置換が行われるのに注意。

で、psファイルまで変換が終わると次の出力が得られるはず。

今15分くらいでパパっと作ったので若干フォントのサイズや種類のバランスが悪いけど、とりあえず置換はできた。簡単！

Rの出力に比べてずいぶんきれいだし、あと数式部分のソースも読みやすく書きやすいと思う。 LaTeXだとドイツ文字のIとか使えていかにもそれっぽいね！

1. すべてのx∈Aにたいしてx≦a
2. すべてのx∈Aにたいしてx≦u ならばa≦u

というようなaは，uをxにいくらでも近づけることができるために，x≦aというようなaの中で最小のもの， Aのすべてを抑えるぎりぎりのものとなる．このとき，これを∪Aと書く．上からギュッと押えるイメージ．これをAの上限と呼ぶ．

先の1.，2.において不等号を逆にしたもの，つまりAを下から押えるぎりぎりのものを∩Aと書き，Aの下限と呼ぶ．

特にQの場合は上限をsupA， 下限をinfAなどとも書く．

上限，下限は必ずしも存在するとは限らない．上限・下限が存在しない例として次のようなものが挙げられる．

• ∪Q　　…　　有理数全体を「抑える」 ことはできない
• ∪φ　　…　　抑える対象が存在しない
• A={x∈Q|x2≦2}    …　　√2は有理数の中に存在していない(Aの中と外の間には途切れ，gapが存在している)

また，P(E)の場合， A∈P(E)という部分集合の集まりについて∪Aは全てのX∈Aを含むもの， つまり合併集合となる． そして∩Aというのは全てのX∈Aに含まれるもの， つまり共通集合となる．また、 Eの部分集合をひとつも取らないΦ⊆P(E)については普通次のようにする．

• ∪Φ=φ
• ∩Φ=E

これは最初の1.，2.に当てはめてみると，1.についてはそもそも判断すべきXが存在しないので自動的に満たされ、 2.についてはUというのが任意のU∈P(E)になるため、「任意のUに含まれるX=φ」 「任意のUを含むX＝E」というようなことになる。ただ、なじめなかったら定義と思ってしまってもいいそうだ。

また，Q，P(E)について次の性質が挙げられる．

• ∪{x|x<a} = a
• ∪{x|x>a} = a

ただし，自然数の集合Nでは次のようになって成立しない．

• ∪{x∈N|x<n+1} = n

関数 f:Λ → Q について， たとえば∪f(Λ)，∩f(Λ)というようなものを表したいとき，次のような書き方が使われる．

一般のxnについては同様に次の記法が使われる．

（ちなみにLaTeXで∪のでかいやつを描くコマンドが\bigcupだと調べるのにちょっと苦労した． 考えたらわかりそうな気もしたけど．）

また，A⊆QについてAをQに埋め込む単射 i:A → Qではxとi(x)を同一視して次のような書き方がされる．

また，増大列(x0 ≦ x1 ≦ x2...)について，

と

は同値である．

なぜなら，b<a<cについて，n≧Nならb≦xn≦cとなるようなNがとれる．いま考えているのは増大列なので， xn≦cはすべてのnで成立する．そこでcの下限を考えると， xn≦aとなる．

----

まだ大丈夫だけどこれは他の本(高木貞二の解析概論とか．ほんの1mmくらいかじっただけだけど．)で見たからだな． どうもやっぱある程度高等な予備知識を要求されている気がする．

このままだとすぐ付いていけなくなりそうなのでなんか別にもう一冊くらい探さないと．などと思って集合・ 位相入門 を頼んだ．ちょっと頑張ろう。

教科書

• 森 毅　『位相のこころ』 (ちくま学芸文庫)

順序構造を考える対象

（O1） x ≦ x

（O2） x ≦ y,　y ≦ z　なら　x ≦ z

というような「公理」が成立するもの。Qでいうところの大小関係、 P(E)でいうところの包含関係（≦の代わりに⊆を使うアレ）など。

・「公理」を満たす極端な例

x ≦ y なら x ＝ y　… 　自分以外と比較できない。比較ができるってことは、そいつが自分っていうこと。 離散順序とでも呼ぼう。

任意のx, yで x ≦ y　…　どいつとでも比較していいけど、 どっちがでかい（あるいは同じ）ことにしてもいい。 混沌順序とでも呼ぼう。

こういう変なのは避けたい…

（O1）、（O2）だけでは不十分（擬順序という）。「順序」 というときは次を加える。

（O3) x ≦ y, y ≦ x　なら　x = y

比較のときにxとyを入れ替えてもいいなら、そいつは「同じ」なんだ、という縛りを入れてやる。 別のもの同士を比較するときは入れ替えちゃだめというわけで、混沌順序を避けられる。

また、昔はこれも加えていた

(TO) 任意の x, y は　x ≦ y　または　y ≦ x

けどこれはなくてもいい。 あるとQは大丈夫だけどP(E)とかが順序の仲間に入れない （二つの有理数はいつでも大小関係を定義できるけど、集合XとYはいつも包含関係にあるとは限らない）。 これを満たすものは全順序と呼んでもうちょっとイイモノだとしておく。

次に離散順序を避けるには…

(DO) 任意の x, y について u ≦ x, y ≦v となる u, v が存在する。

集合から適当に2つの元をとってきたとき、 その両方より小さいuとその両方より大きいvが存在するということ。 xとyが別物だったとき、 その両方と比較できるu, vが存在するということだから、 離散順序みたいに比較できる集まりがポツポツと飛んでいるものは避けられるという仕組み（と思う）。

これは有向順序という。Qみたいな全順序はもちろん有向順序。 P(E)もその元Xについてφ⊆X⊆E（少なくとも空集合は含むし、空間には含まれる） なのでやっぱり有向順序。

ついでにもうちょっと厳しくする。

（LO) 任意の x, y について

a ≦ x, y ≦ b

u ≦ x, y ≦ v　なら　u ≦ a ≦ b ≦ v

となる a, b が存在する。

a, b は次のように書く。

a = x ∩ y

b = x ∪ y

具体的に言うと、

aはxとyの両方より小さい、または両方に含まれる、 といったようなもののなかで一番デカイやつ

bはxとyの両方より大きい、または両方を含む、 といったようなもののなかで一番ちっこいやつ

aは下から押さえる感じ（だから、「デカイやつ」とか言っても実態は普通小さい。小さい中で一番デカイ）。bは上から押さえる感じ （同じく「ちっこい」といっても実態はでかい）。こういうものを束順序という （ちなみにこの束(laticce)はファイバー束とかの「束(bundle)」とは違う。違う概念に同じ名前。ややこしい！）。

ただの言い直しだけど、

Qの中の話だったらaはxとyのうちで小さいほう、bは大きいほうのこと。

P(E)の中の話だったらaは集合XとYの共通集合（共通部分）、bは合併集合（和集合）。

…というように、「順序」 と一言に言っても向きもへったくれもないものや自分と違うものとは比較を拒むようなものからいろんな段階のものが考えられるんだけど、 適当にx, ｙと2個持ってきたときに、 それらをまとめて押さえられる何かを定義できるようにしておくと便利ですよと、そんな話（だろうと思う） 。それらをまとめて押さえられる何かを定義するのは次のパートですね。

----

位相とか集合とか一度ちゃんとやらないと、とか張り切って読み始めた。

…

1時間で3ページしか進んでないとかどういうことなの…。

a ≦ x, y ≦ b　とか書いてあるときに「a ≦ x」, 「y ≦ b」と勘違いしてて30分くらいロス。 xとyまとめてるだけなのに。x ≦ y,　y ≦ zと定規でもって比べたら、a ≦ x, y ≦ bの方はカンマの後のスペースが0.5ミリくらい短かったよ！

本の最初のほうはその人の記号の使い方とか説明の仕方とかに慣れてないからどうでもいいとこで躓いてしまう。 この辺の話だって他の本で見たはずの内容だけどあらためて丁寧に追ってみると今みたいなつまらない躓きが…。 まあでも慣れたからペース上がるはず。きっと。おそらく。多分…。

ところで（TO）とか（LO）とかいう書き方は何なんだろう。人名とかの略かな？

教科書

• 森 毅　『位相のこころ』 (ちくま学芸文庫)

要するに、メンデルの法則が正しいかどうかはさておいて、メンデルの実験には捏造疑惑があるのだ（もちろん、 メンデルの法則がいくつかの例外を除いて、 あるいは例外についても適当な説明方法によって正しいということが主張できるのは周知のとおりだが）。

しかし、日本植物整理学会の「みんなのひろば」という一般の人からの質問に答えるコーナーの質問と回答をまとめた本、 「これでナットク! 植物の謎 」 という本で、次のような記述があった。

しかしメンデルにならって交雑の際に慎重に個体を選んだ三人の研究者が、それぞれ独自に出したデータも、メンデルと同様に 「きれい」なものだったので、フィッシャーの批判はあたらないとされています。

本当か？

ちなみに、Webページはこちら。 編集の都合だろうが書籍よりも記述がいくらか多い。

おそらく、「三人の研究者」というのはメンデルの法則を再発見したド・フリース、コレンス、 チェルマクの3人のことだろう。「工学のためのデータサイエンス入門」 にメンデル、コレンス、チェルマクのデータが載っているので引用する。

実験者　　黄色個体　　緑色個体
メンデル　　　　6022　　　　　2001
コレンス　　　　1394　　　　　 453
チェルマク　　 3580　　　　　 1190

ここでそれぞれのp値は0.9076, 0.6480, 0.9467である。

p値というのは、「理論が正しいとして、この比率よりも外れた比率が観測される確率」を意味する。つまり、 たとえばメンデルのデータから観測された0.9076という値は、「理論が正しい」ときに観測されたとしても問題のない範囲に存在している。

と、ぱっと見はそう思うだろう。

ただしこれを逆に考えると、「理論が正しいとして、この比率よりも『理論に適合する』 比率が観測される確率は1-pである」ということが言える。つまり、「理論」、3：1という比率、 これにメンデルの観測データはきわめて近い。そしてメンデルのデータよりもより「正確な」 データが得られる確率は1-0.9076=0.0924にすぎない。

コレンスはともかく、チェルマクのデータについてはメンデルよりも「正確」であり、 このデータが得られる確率は1-0.9467=0.0533、つまり5%程度である。

さらに悪いことに、メンデルが予想した「誤った理論」に「正確」なデータもメンデルは得てしまっている。 そこではある対立遺伝子についてF1個体の分離比をAA:Aa＝201：399と報告していたのだが、 AAとAaの個体の識別方法に問題があった。メンデルは優勢形質が得られた個体からの種子10粒を育てて自家受精し、得られた種子にaa型 （劣勢）形質が見られれば親はAa型であり、全てAA型（優勢）ならば親はAA型であると判断を下した。 しかし確率的には(3/4)^10≒5%で親がAaにもかかわらず子が全てAAになってしまうという場合が生ずる。 よってメンデルの方法によって観測されるべきF1世代の遺伝子型の分離比は223：377でなければならないのである。

それに加えてなお悪いことに、メンデルは他にも数多くの形質について分離比を調査しているのである。そして、 フィッシャーの批判というのはこれらのデータ全てを総合して考察した結果にたいして向けられている。フィッシャーの計算するところによれば、 メンデルが報告した値よりも「正確な」データが得られる確率は、0.00005、 つまり10万回に5回という極めて小さな値となるというものだった（これについては「やさしい統計入門 」 から引用）。

以上から次のことが言える。

まず、メンデルのデータが当てはまりすぎているのは事実である。

つぎに、メンデルに続いた「三人の研究者」のデータが「当てはまりすぎている」かどうかというと、 これは微妙なところだ。ド・フリースのデータについては調べていないので分からないが、とりあえずチェルマクのデータはたしかに「きれい」 に当てはまっている。しかし、コレンスのデータはp=0.6480であり、これは実際の実験で「ありそうな」値だ。決して 「当てはまりすぎている」というものではない。

なので、続く研究者のデータも「きれい」なものだったからフィッシャーの批判はあてはまらない、 というようなことは言えない。

もちろん、植物が確率に従わないような何らかの制御機構を持っており、分離比が確率の示す以上に3： 1に近づくのだ、というような説は否定できない。ただ、「工学のためのデータサイエンス入門」には 「三人の研究者」 以外の研究データもいくつか載っており、 それらのp値だけを挙げると0.7408, 0.3779, 0.1742, 0.4488, 0.0950である。 理論に対して適当なずれを持っており、「きれい」というわけではない。

どうも私の調べた範囲ではフィッシャーに分があるように思えるがどうだろうか。

ちなみに、フィッシャーの意見は 「どのような結果が得られることが望ましいかを熟知している数名の研究補佐員によって、メンデルは欺かれていた可能性がある」というもので、 メンデルを直接に攻撃するものではなかったようだ。

つ「っぽい写真」

受付にいた。

イチゴ狩りとか生まれて初めてだ。もうこんな機会はないかもなぁ。

今回は人間が多く移っていたのでFlickrに上げるほどのアレはない（あまり面白くない写真を数枚上げたけど）。 それになんか結構楽しかったので、今後語る機会はない。

ただ、デジカメのオートブラケット機能は常時ONにしておくべきだということはよく分かった。 8～16GBのメディアが数千円で買える今、容量による制限はない。無駄に複数の露出で無駄に連写しておくと、「たまたまうまくいった」 写真が結構入る。一枚手ぶれしても他の数枚にブレなしの写真が記録される可能性は高い。 シャッタータイミングとかもあまり気を使う必要がない。買ったときからONにしとけばよかったと悔やまれる。

定義の仕方はいろいろある。けども、大体以下のポイントを抑えたものが使われる（っぽい）。

• -∞～+∞の区間で積分すると1になる
• x=0において無限大の値をとる

前回の(1)の定義によるとδ(x)はx=0以外の点で0にならないといけないんだけど、これはどうも無視されることがあるらしい。 というのも、デルタ関数を実際に使うにあたって重要なのは前回でいう性質の2つ目、式(4)のところで説明した「値を抜き出す」 ということであって、これが守られれば少々のことは無視してもおｋということらしい。

とりあえず教科書に載ってるのを列挙する。過程とかあまり書いてないので説明は無理っぽ。 ただ少なくとも最初に書いた2つの点は抑えてる。

1個目。

2個目。

これは実質的に1個目と同じらしい。Imというのは複素数の虚部を取り出す演算子。iは虚数単位で√-1。

3個目。

おなじみ正規分布の密度関数を平均=0のもとで分散を0に限りなく近づけた場合と同じ。 これはもともと確率の密度を与える式なんだから全区間積分すると1になるってのは分かりやすい。 分散が0になるとx=0の部分が無限大に近づくというのもイメージしやすい。個人的に一番「それっぽい」気がする。

4個目。

これが最初に言った定義を「少々無視」する形式。sin axはaをでかくすると0付近でどんどん激しく振動するだけで収束はしない。 ただ、この形をデルタ関数としておくとフーリエ変換の際などに便利らしく、この定義が結構頻繁に使われるぽい。

教科書的には式(4)は三角級数からいくつかのステップを踏んで誘導する感じで紹介されてるけど、 そもそもその三角級数が出てきた経緯がちょっと分からなかったので完成したブツだけ。

この4つについて適当なaの値を与えてグラフをプロットし挙動を観察してみる。

まず式(1)。そんな無茶苦茶早く収束するっていうわけでもない感じ。

式(2)。プロットすると式(1)と同じものだってことがよく分かる。ちなみにRで複素数の虚部を取り出す関数はIm()。 複素数のオブジェクトを作って虚部にだけ値を与えるにはcomplex(imaginary=1)などとする。 一応ここだけプロットに使った関数のソース。

Delta.2 <- function(a, x){
 a <- complex(imaginary=1/a)
 Im((pi*(x-a))^-1)
 }

式(3)。収束が早いっぽい。

そして式(4)。aの値を大きくするとブレがむしろ拡大していくのが分かる。

※「ブレが拡大」という部分ですが、軽い指摘を頂いたのでよくよく考えてみたら確かに微妙な表現だと思えてきました。 x=0の付近で値は激しく振動し、ピークは鋭くなっていくわけですが、周辺での振幅は拡大されるわけではありません。 試しにa=200とした場合のグラフを描いてみました（↓に追加してあります。 もっと大きくてもいいのですがこれ以上だと塗りつぶされてしまいます。）。周辺部の振幅で言えばこれは同じと言うべきでしょう。重要なのは 「収束はしないものの、激しく振動しつつx=0でのピークが高く鋭くなっていく」という点でしょう。ちなみに教科書↓の表現を要約しますと、 「sin axはa→∞としても0となるわけではない。しかしa→∞で激しく振動するため、『値を抜き出す』性質は成立する」。 (2009/02/24 追記)

で今日はおしまい。ええ、手抜きです。

…なんだけど、Rで数式書くのが案外面倒だった。なので数式の部分だけメモっておく。 使うときはプロットのmain=とかの後にそのまま突っ込むだけ。

#式(1)
 expression(lim(over(a, pi(x^2+a^2)), a %->% 0))
#式(2)
 expression(lim(paste(Im, over(1, pi(x-ia))), a %->% 0))
#式(3)
 expression(lim(paste(over(1, paste(a, sqrt(pi))),exp,
 bgroup("(", paste("-",over(x^2, a^2)), ")")), a %->% 0))
#式(4)
 expression(lim(over(paste(sin," ",ax), paste(pi,x)), a%->%infinity))
#ついでにプロット時の余白の設定
 par(mar=c(3,3,3,1), mgp=c(2,1,0))

試行錯誤しながら書いた部分があるので正しくないかもしれない。式(3)とかなんだかpaste()たくさん使ってるし。 expression()の使い方あまり把握しきれてない気がする。

やっぱ数式書くならLaTeXのが楽だ。いちいち書き方調べなくても勘で大体書けるし。 LaTeX使って式を挿入する方法を調べないといけない。 どうもPSfragとかいうのでポストスクリプトの中身を書き換えられるらしい。 それでいけるのかも。また今度。

教科書

• 小野寺嘉孝 『物理のための応用数学』 （3.3 デルタ関数のいろいろな表現 p.50～ ）

基本的に0だけど、x=0の点を含む積分では1になる。積分範囲は0から任意の実数ε （0とほとんど同じくらい小さくても、無限大でもいい）だけ離れた範囲。x=0の点で無限大の値を持つ、 というような表現がされることもある。

本来積分を含めてしか性質が現れてこなくて、そういう意味で「超関数」とか呼ばれてる（多分）のだけれど、 無理やりグラフにするとたとえばこんな感じ 。

εはすっごい小さい数で、とりあえず形は正規分布になってる。正規分布において平均=0、 分散→0の極限をデルタ関数として定義することもできるらしい。

まあとにかく、デルタ関数は普通0を中心に対称な、偶関数として定義される。よって、 偶関数として定義されたなら、という条件付で

とりあえずここまで定義。以下特徴の説明。デルタ関数が導入されるときによく説明される「性質」をいくつか。

それで、今の図において－ε～εの範囲で線の下の面積が1なんだから、－ε～0、 0～εの範囲に含まれる面積は1/2じゃないといけないだろう。よって、定義から導かれる性質の1つ目。

これ以上特に説明も不要だろうから次へ。

性質の2つ目。

任意の関数f(x)とδ(x-a)（この場合デルタ関数はx=aにおいて値を持つ）の積を積分する（区間はx=aを含む）と、 x=aにおけるf(x)の値が得られる、ということを言っている。 これは次のように説明できるだろうと思う（表現として正確じゃないかも）。

積分範囲をx=aから±εの範囲とする。εは任意の数なので、非常に狭い区間を考えるとf(x)→f(a)。 右辺δ(0)×dxという部分は、リーマン和の一切れだとイメージしてほしい。 縦×横で面積が出てくるということを言っている。縦がδ(0)、横がdx。それで、デルタ関数の下の部分の面積というのは1だったから、 結局f(a)が残る。

次、性質の3つ目。

「中」をa倍するのは、「外」から1/|a|倍したのと同じこと、ということを言っている。 これは次のように証明する。

まず、δ(y)というものを考えると、定義より∫δ(y) dy = 1。ここで積分変数をyからaxへ変更する。置換積分だ。

dy/dx = aなのでdy=a×dx。よって次式が得られる。

2段目の＝の後は変換ではなくて単に等しいということを言っているだけ。定義から積分の結果はどっちも1になるので。それで、 式(7)の後半2つを使うと次が得られる。

一段目は式(7)をaで割っただけ。2段目は両辺をdxで微分している。積分して微分したら元に戻るので単純に積分記号をはずすだけ。

しかしaを-aとおいて式(8)まで誘導すると1/aが-1/aとなってしまう。

ここで最初の定義(2)を思い出すと、偶関数なのでδ(ax) = δ(-ax)。よって、aを-aとして計算した場合では、 式(7)の時点で得られる -a∫δ(-ax)dxの積分の中身が必ず正となるため、 a < 0である必要性が生じる。逆に言えば、 a < 0であるとき、 式(7)の時点でaの前に-をつけなければならない。

よって-aからはじめて式(8)の時点で得られる-1/aという値は非負でなければならない。

なので式(6)のように絶対値記号のついた表現となっている。

次。性質その4。

性質3の拡張なのだけれど、axをいきなりf(x)としているから分かりにくい。aiという定数が右辺に含まれるが、 これは関数f(x)の実根である。つまり、f(ai) = 0となるような値。

ちなみにf(x)は任意の関数というわけには行かない。右辺でゼロ割りが生ずる危険性があるので、 f(x)はf'(a)≠0という制限がつく。これは結構厳しい制限だ。 x軸と平行に交わる関数はすべてだめなのだから。たとえばx^2という単純な関数を考えるとその根は0、 導関数は2xなので、f'(a) = f'(0) = 2*0 = 0となってしまうために式(9)は適用できない。

さて、式(9)は根が複数あるような関数を想定しているが、とりあえずひとつからはじめよう。

まず式(10)は部分積分。yをf(x)とおいて積分変数をxへ変更している。

次、この積分が非0となるのはデルタ関数の中身が0となる点を含む区間においてのみであるから、 f(x)が0となるようなx、 つまり実根aから±εに積分範囲を限る。この操作を表しているのが式(11)。

εは任意の実数なのでいくらでも小さくとることができる。εを十分小さくとったとき、 f'(a±ε)をf'(a)という定数で近似できる。 定数は積分記号の外へ括りだせる、 ということを言っているのが式(12)。

式(13)は式(7)のときと同じで変換ではなく、単に等しいということを言っている。式(11)で積分範囲をa±εに限ったので、 xがそのような区間しか動けなかったとしてもちゃんと値が1になるようにしなければならない。 そのためにδ(x)ではなくδ(x-a)としている。式(13)の積分範囲はa±εより広ければ無限大でもいい。 ちょっとずれてるだけで基本のデルタ関数だ。

それで、式(12)と式(13)から次が得られる。

そして、実根が複数ある場合を考えると、 どの実根を含む範囲においても式(14)が成立する必要がある。 これを表現するにはすべての場合を単純に足し合わせればよく、したがって次式が得られる。

そして、式(6)で絶対値記号をつけたのと同様の議論をすることにより、やはり絶対値記号つきの式(9)が得られる。

----

教科書なんか見てると四つ目の性質、もしくはそれに類似する「性質」や「特徴」の説明って結構証明なしに与えられることが多い。 もしくは証明することが課題になってたりして。しかしΣと絶対値記号が突然出てくるのは結構ビビる。混乱する。

そんなのちょっと考えればすぐに解決できる程度のレベルの人間を想定してるのか、 あるいはいい先生に教えてもらうのが普通なのかは分からない。

でも式(9)の形の変換ってグリーン関数とか使う場合に出てくる気がする （伝熱関係の論文で見たことある気がするだけで使ったことはないんだけど）ので、ちゃんと把握している必要がある事項だと思う。

ここに躓いたのいつだっけか。確か半年くらい前。そのときは何が道具として必要なのかの見当がつけられずにあきらめたけど、 結局は置換積分と「εをすっごい小さくする！」という微積分によくある操作だった。

やっぱ基礎だな。基礎を徹底的に固めないといけない。誰か固めるべき基礎をリストアップしてくれ。

----

参考文献

• 小野寺嘉孝 『物理のための応用数学』 （3．デルタ関数 p.47～ ）

http://movapic.com/nozma/pic/130840

味付けとしてはだしの素としょう油だけ。

ご飯を水でよく煮て、刻んだニラを加えて、さらによく煮て、卵入れて完成。

風邪引いたときには必ず出てきたけどニラが常備されていた覚えはない不思議。

作り方というほどの作り方もないけど、一応ばーちゃんに問い合わせて確認したので忘れないようメモ。

確か家ではニラのクロロフィルが変性して若干色が変わるくらいまで煮ていた気がする。味は時々薄かったし、 ニラの量もまちまちだったので分量らしい分量はないはず。だしの素としょう油だけなので薄かったらしょう油かけたら大丈夫。

とりあえず100円で買ってきたニラの束が1/3も減ってないのでしばらくニラおじや食う。

というのが最小作用の原理。もしくはハミルトンの原理。積分Sは作用と呼ばれる。Lは一般化座標q, その導関数q'（一般化速度） , 時間tからなる系を特徴付ける関数で、ラグランジアンという。（「特徴付ける」って何なのさ。という感じだけど、とりあえずここは 「ラグランジアンとして具体的な形を考えなくてもしばらく計算が進められる」 ということを強調するためにこんな微妙な言い方をしているだけなので、細かいことは考えず先へまいりましょう。）

一般化座標というのは系において位置を決めるのに必要な任意の個数の量の組で、たとえば平面上で垂直方向へy、 水平方向へxと座標を設定すれば、一般化座標の一例は(x, y)となる。「一例」 といったのはたとえば水平面からの仰角θと直線距離rを用いて(θ, r)としても位置は一意に決定できるからだ。 ちなみに位置を一意に決定するのに必要な量の個数をその系の自由度と呼ぶ。

さて、積分Sが最小となるということの意味を数式で表してみよう。変分記号δを使えば、それは次のようになる。

ここから前に変分法を説明したときと同様の手順を踏んでオイラー方程式を得る。

このオイラー方程式はラグランジュ方程式と呼ばれる。

次にラグランジアンの具体的な形を示す。重力のある質点系においてラグランジアンは次の形で定義される。

Tは運動ポテンシャル。Uは位置ポテンシャルを表す。

いきなり出てきた感があるけど、これは相互作用のない系におけるラグランジアンに、相互作用（質点ともうひとつの質点、 地球との相互作用である）によって決まる座標の関数Uを足し合わせることによって得られている。

運動ポテンシャルは力（F = ma = m dv/dt)を位置について積分することで得られる。 物理の教科書とかによく載っているので多少飛ばし気味で（ただしランダウ=リフシッツ「力学・場の理論」では以下と異なり、 先に述べたように「相互作用のない系におけるラグランジアン」として運動ポテンシャルを導いており、議論は最小作用の原理だけを用いて進む。 が、ちょっとまだ理解しきれてない部分があるのでそれはパス。気になる人は「力学・場の理論」でどうぞ！）。

一列目から二列目の変換は置換積分であることに注意。 ここで積分変数がqからvへ変わっている。

位置ポテンシャルも同じノリで重力（F = mg)を高さhについて積分することで得られる。

こっちは置換積分とかしない。

式(5)、(6)から重力のある質点系での具体的なラグランジアンの形を求めれば次の形になる。

これをラグランジュ方程式(3)へ代入すると、次の式が得られる。

ところでここで左辺第一項はma（aは加速度)である。第二項、計算前の形で∂U/∂qであるこれは力である。

そう、要するにこれはma = Fというニュートンの運動方程式に他ならない （別にニュートンの運動方程式を求めるだけならUとして具体的に位置エネルギーというものを挙げず、 座標と相互作用により決まるポテンシャルエネルギーとしておけばよかったんだけど抽象的すぎて分かりにくいので）。

たとえばここで垂直方向へy、水平方向へxとして座標を設定し、式(8)へ座標の情報を取り込むと次の式が得られる（空気抵抗は無視） 。

これはよくある（と思う）空気抵抗を無視した状態で物体を投射したときの運動を求める問題。まあ放物線が出てくるわけですが、 結構面倒なんで今日は解かずにここで終わり。

とりあえず、最小作用の原理からニュートンの運動方程式が出てくる、というお話でした。

----

参考文献

• 力学・場の理論―ランダウ=リフシッツ物理学小教程 (ちくま学芸文庫) （最小作用の原理だけからニュートンの運動方程式を導く方法などについて。第1章。）
• 新・物理入門 (駿台受験シリーズ) （放物線を求める問題についてはp.25）