2008年01月30日

ベクトルの水準に関する覚書[R]

問題は存外むつかしいようだから、出鱈目に試行して失敗を繰り返していては、なかなか解けそうにない。 敵は思いのほかに手剛いから、作戦計画を練って、取り掛からないと、辛い目を見るであろう。おっくうだけれども、 少し頭を働かせなくてはならない!

Rにおいて文字列のみを要素とするベクトル、例えば

x <- c("A0", "A0", "B1", "B1", "C2", "C2")

に対しては、関数factor()を用いることで要素をグループ化することができる。

x <- factor(x)

とやると

> x
[1] A0 A0 B1 B1 C2 C2
Levels: A0 B1 C2

このように、3つの水準が定義される。だが、ベクトルを作るときに間違えて

x <- c("A0", "A0", "B1", "B1", "C2", "C2")

とやってしまったとしよう(2つ目のA0の0が全角になっている。

その後

x <- factor(x)

とやると

> x
[1] A0  A0 B1  B1  C2  C2
Levels: A0 A0 B1 C2

もちろんこうなる。水準が一つ余計だ。というわけで、修正を加える。

x[2] <- "A0"

ベクトルxの2つ目に"A0"を入れた。これで水準は3つになるはず。が、

> x
[1] A0 A0 B1 B1 C2 C2
Levels: A0 A0 B1 C2

既に無い水準が定義されたままになってしまっている。

この些細なミスがときに深刻な影響をもたらすことがある。分散分析や多重比較など、 水準を用いる関数ではとくに注意が必要。

対策は簡単。ベクトルに修正を加えたら

x <- factor(x)

を再度実行すればいい。

昨夜はこの単純なミスに気付かないばかりにかなりの時間を費やしてしまった。 分散分析も多重比較もプロットも出来ないのだから相当あせった。

ところで、分散分析をしてから多重比較というのはダメなのだろうか。

いや、正しくないには確かに正しくない。分散分析はF統計量を使うのだし、 それが多重比較のステップに組み込まれていないのであれば多重性の問題が生じてしまうという主張(永田・ 吉田「統計的多重比較法の基礎」) はもっともだ。でも、多くの学術誌が段階を踏めと言っているし、郡分け変数の影響を知りたいにしろ、 郡間の差を知りたいにしろ、 どこかの郡にはっきりとした差があれば結果は変わらないのだから当面は段階を踏んでおけばいいんじゃないのか(中澤港 「Rによる統計解析の基礎」)と言われればなるほどなと思う。

大体、多重性の問題がどのように発生して、具体的に危険率が何%になってしまうのか検討が付かない。 t検定を繰り返したときの危険率もいまだに良く分からない。例えば3標本にt検定を繰り返す場合、 データが独立ではないので3回繰り返しても単純に1-0.95^3とはならないみたいだ。 実際は1-0.95^3よりも少し少ない程度に収まる。

分散分析と多重比較だってそれぞれの危険率はある水準以下に抑えられているんだから、 多重性の問題が発生したって大した問題にはならないんじゃないだろうか。まあ、大した問題でないから無視するというのもアバウトな話だが。

というか、自然科学に関するデータならば最終的には「見て」判断しなければいけない。それに、 棄却検定である以上データ数を増やせば処理間に差は必ず出る。この世の中に対象に微塵の影響も及ぼさない処理などあろうはずもない。 自然を相手にする以上もとからしてアバウトなのだ。考慮するなら分散分析と多重比較の間の多重性よりも、 サンプルサイズと検出力の関係だとか、そもそも棄却検定を適用することについての妥当性だとかそのへんが先じゃないだろうか。 こういうことをフィッシャー先生はどこまで考えていたんだろう。図書館に本あったし借りてくるか。

とか昨日人の話を聞いていて思った次第。

posted by Rion778 at 02:05 | Comment(0) | TrackBack(0) | PC関連。HTMLとか,Linuxとか,Rとか | このブログの読者になる | 更新情報をチェックする
2008年01月26日

単位の異なる複数の時系列データを一つの表にプロットする[R]

どのデータも同一の時系列に沿って採取しているのだが単位が違う. そういう時はグラフの右と左で違う目盛りを使ったりして対応するだろう.

しかし,データが3種類を超える場合はどうしようもないので別々に作ったグラフを縦に並べて比較を試みるだろう.例えばこんな風に.

3graph4

今日 はこのグラフをRで作ることを目標とする.

まずはダミーのデータを用意.

days <- 1:30
dry  <- sqrt((rnorm(10))^2)+sqrt(days)
fresh <- dry^2+10
ratio <- dry/fresh

まーダミーなので解釈は何でもいいんだけど,一応daysがある植物を移植した後の経過日数,dryが乾物重,freshが新鮮重, ratioが乾物率(乾物重/新鮮重)という設定で.1ヶ月も毎日サンプリングして新鮮重と乾物重なんか測定してたら気が狂いそうだな….

まーとにかくまずは普通に描画領域を3分割してプロットしてみよう.

par(mfrow=c(3,1))
plot(days, fresh, type="b",ylab="Fresh weight(g)",
 xlab="Day after transplant")
plot(days, dry, type="b", ylab="Dry weight(g)",
 xlab="Day after transplant")
plot(days, ratio, type="b", ylab="Dry weight ratio",
 xlab="Day after transplant")

最初のpar(mfrow=c(3,1))で描画領域を縦方向に3分割している. 例えばpar(mfrow=c(3,2))とやったら縦に3分割,横に2分割とかもできる.で, そのあとプロットを実行すると左上から埋まっていく.この場合は上から.何はともあれ実物を.

3graph1

まーこれで も縦にスペースを多くとれば見えないことは無い.時にはこの形式が適することもあるだろう.が, 紙面の節約のためにはグラフ同士をぴったりくっつけたいところ.

そのためには,par(mar=...)というグラフィックスパラメータによってグラフ周りのマージンを指定してやる. mar=の後には4つの数値からなるベクトルを代入する.c(下,左,上,右)という対応なので, 例えばpar(mar=c(1,3,1,3))と指定すれば左右のマージン3行,上下のマージン1行の状態で描画される. ただし左と下は4行以上のマージンが無いと文字がはみ出るので注意.

それで,これを組み合わせて,一番上のグラフは下マージンなし,真ん中は上下マージンなしとして指定してやると

par(mfrow=c(3,1))
par(mar=c(0,4,4,4)
plot(days, fresh, type="b",ylab="Fresh weight(g)")
par(mar=c(0,4,0,4))
plot(days, dry, type="b", ylab="Dry weight(g)")
par(mar=c(4,4,0,4))
plot(days, ratio, type="b", ylab="Dry weight ratio")

3graph2

くっついた.でもくっついただけ.真ん中のグラフだけ大きいし,x軸目盛りは全てのグラフに入っているしでどうにも不恰好.

まずはx軸の目盛りを消す.それには,プロットするときに引数としてaxes=Fを指定する.そうすると軸目盛りは全部消える. それどころか周りの枠線も消える.そのままでは困るので,その後すぐにaxis(2)と打ち込む.するとy軸の目盛りが書き込まれる. ついでにbox()と打ち込むと周りの枠も復活する.

グラフの個数が2つならばこれまでの方法で2つのグラフを1つにまとめても問題は無いだろう.

しかし,3つ以上つなげたときに高さが違ってしまうのは防げない.

ので,全てのグラフの上下マージンを消してしまう.つまりpar(mar=c(0,4,0,4))と最初に指定したら, 全てその余白設定でグラフをプロットする.

そうすると上下マージンがなくなってしまうが,外周マージンを指定してやれば問題は解決する. 外周マージンはpar(oma=...)と指定する....の部分は先のマージンの時と同じく,ベクトルとして4つの値を入れてやる. 上下に4行分のマージンが欲しければpar(oma=c(4,0,4,0))のように.

しかしこれだとx軸のラベル名やタイトルが領域外となってしまうので表示されなくなる.そこで,mtext関数を使う.例えばx軸に 「Day after transplant」というラベル名が欲しい場合,

mtext("Day after transplant", side=1, line=3, cex=0.8)

とやると,普通にラベル名を指定したのと同様の場所にテキストが表示されるはずだ.sideはテキストを表示する位置(1:下,2: 左,3:上,4:右)の指定,lineは領域から何行離すか(3行空けでラベルと同様)の指定, cexは文字サイズの指定(値は標準に対する倍率).なので,タイトルが欲しい場合は一つ目のグラフを描画した直後にside=3, cex=1くらいの指定でmtext関数を実行するとそれっぽくなると思う.

まとめよう.

1. 複数のグラフを同時に表示するには領域を分割し,その後プロットする.

par(mfrow=c(3,1))
plot(...)
plot(...)
plot(...)

2. 複数のグラフをぴったりくっつけるには上下のマージンを0にする.全体としてのマージンは外周マージンとして指定する.

par(mfrow=c(3,1))
par(mar=c(0,4,0,4))
par(oma=c(4,0,4,0))
plot(...)
plot(...)
plot(...)

3. x軸の目盛りを消すには,プロットの引数としてaxes=Fを指定したのち,axis(2)とbox()を実行.

par(mfrow=c(3,1))
par(mar=c(0,4,0,4))
par(oma=c(4,0,4,0))
plot(... , axes=F)
axis(2)
box()
plot(... , axes=F)
axis(2)
box()
plot(...)

4. x軸ラベルはmtext関数で入れる

par(mfrow=c(3,1))
par(mar=c(0,4,0,4))
par(oma=c(4,0,4,0))
plot(... , axes=F)
axis(2)
box()
plot(... , axes=F)
axis(2)
box()
plot(...)
mtext("...", side=1, line=3, cex=0.8)

以上のポイントを押さえれば冒頭に示したようなグラフが描けるはず.

ちなみに冒頭のグラフを書くのに使ったスクリプトはこちら.

par(mfrow=c(3,1))
par(mar=c(0,4,0,4))
par(oma=c(4,0,4,0))
plot(days, fresh, type="b",ylab="Fresh weight(g)", axes=F)
box()
axis(2)
plot(days, dry, type="b", ylab="Dry weight(g)", axes=F)
box()
axis(2)
plot(days, ratio, type="b", ylab="Dry weight ratio")
mtext("Day after transplant", side=1, line=3, cex=0.8)

でも時系列データを取り扱うパッケージで同じようなグラフを見たことあるような気がするので, もしかするともっと楽に書く方法があるのかもしれない.

…ってグラフィックデバイスの使い方をメモるのを忘れてた.

pngへ出力するには,最低限

png("ファイル名.png", width=... , height=... )

の指定があればいい.widthとheightはピクセルで指定.アウトプット先は現在のディレクトリ (ファイル→ディレクトリの変更で確認).描き終わったら

dev.off()

でデバイスを閉じるのを忘れないこと.

posted by Rion778 at 20:37 | Comment(0) | TrackBack(0) | PC関連。HTMLとか,Linuxとか,Rとか | このブログの読者になる | 更新情報をチェックする

Rでノーマル・チップスを作ってみる

推計学のすすめ― 決定と計画の科学」 という統計学入門書の71ページあたりに,「ノーマル・チップス」 というものを使った遊びが載っている.

ノーマル・チップスというのは0〜60の数字を枚数を変えながら998枚の紙切れに記述することで,近似的に平均が30, 標準偏差が10(分散100)となる集団を作ったもののことだ.これを母集団とみなすことで,標本集団の平均値のばらつき方だとか, 分散比のばらつき方だとかを観察する.

もっとも今ならExcelなんかで簡単に乱数が作れるので, わざわざこんなことをして平均値や分散について講義するようなことも無いだろう.

それでも,実際に数字の書かれているものを混ぜて,取り出して,集計して,というのはイメージがしやすい. PC上のシミュレーションであっても,その感じは少しくらい残るだろう. というわけで暇つぶしがてらRでノーマルチップスを使った遊びを再現してみた.

まず,ノーマルチップスの作成.

norm.tips <- rep(c(0:60),
 c(1,1,1,1,1,2,2,3,4,4,5,7,8,9,11,13,15,17,19,22,24,27,
 29,31,33,35,37,38,39,40,40,40,39,38,37,35,33,31,29,
 27,24,22,19,17,15,13,11,9,8,7,5,4,4,3,2,2,1,1,1,1,1))

…なんてめんどくさいんだ.でもこれで998枚.ちなみに,rep関数は一つ目の引数で数字の種類を, 二つ目の引数で何回繰り返すかを指定する.両方を同じ長さのベクトルで指定してやると, 上のように大量の数値を比較的少ない労力で入力できる.

まずはヒストグラムを描いて,本当に平均30,標準偏差10の正規分布を近時しているのか確認.

library(MASS)
par(mar=c(4,4,2,2))
truehist(norm.tips, ylim=c(0,0.05), ylab="", xlab="",main="ノーマル・チップス")
par(new=T)
curve(dnorm(x,mean=30,sd=10), ylim=c(0,0.05), ylab="相対度数", xlab="チップの番号")

normtips

重ねて描いた曲線が平均30,標準偏差10の正規分布の確率密度関数.たしかにノーマル・チップスは平均30, 標準偏差10の正規分布を近似しているようだ.

それじゃあこのなかから10個をサンプリングする.データが10個たまったらまずは平均値を計算する.そして, 平均値を計算したら再び10個サンプリングし,平均値を計算し…という作業を2000回繰り返す.なお, 10という数字にも2000という数字にもとくに意味は無い.箱の中にノーマル・チップスが入っているとして,よく混ぜて1枚取り出し, 値をメモったら箱に戻し,またよく混ぜて1枚取り出し,値をメモって箱に戻し…という作業を何度も繰り返すことを想像して欲しい. 実際にやったら大変な労力だ.が,Rなら次の5行で実行可能だ.

sample.mean <- numeric(0) #標本の平均を入れる入れ物
for(i in 1:2000){         #2000回の繰り返し
 x <- sample(norm.tips, 10, replace=T) #重複を許して10個サンプリング
 sample.mean[i] <- mean(x)             #標本の平均を計算し,入れ物へ
 }

これで2000個の平均値が"sample.mean"という入れ物に入った.ノーマル・チップスの平均値,分散と, この"2000個の平均値"の平均値と分散をそれぞれ計算して比較してみよう.

> mean(norm.tips)  #ノーマル・チップスの平均
[1] 30
> mean(sample.mean)  #標本の平均の平均
[1] 30.06575
> var(norm.tips) #ノーマル・チップスの分散
[1] 99.18957
> var(sample.mean)    #標本の平均の分散
[1] 9.502943

ノーマル・チップスの平均と,標本の平均の平均がほぼ等しいことはすぐ納得できるだろう.分散の方に注目すると, 標本の平均の分散はノーマル・チップスの分散に比べおよそ10分の1になっている.10というのはサンプリング数に等しい.じゃあ, サンプリング数を5個にしたら5分の1に,20個にしたら20分の1になるのか.なるのだ (上のスクリプトを少しいじるだけなので確認してみてほしい).それどころか,驚くべきことに母集団の分布がどんな形でも, 分散と平均の存在する分布であり,サンプリング数がある程度大きいという条件が満たされれば, サンプリングの平均値は平均値が母集団に等しく,分散は「母集団の分散/サンプリング数」に等しい「正規分布」 をするのである(関連:中心極限定理).

ということは,この「標本の平均の平均」は平均値30,分散10の正規分布をするはずだ.ヒストグラムを描いてみよう.

par(mar=c(4,4,2,2))
truehist(sample.mean, ylim=c(0,0.13),xlab="",ylab="") #標本の平均値のばらつき方をヒストグラムで見る
par(new=T)
curve(dnorm(x,mean=30,sd=sqrt(10)),ylim=c(0,0.13), xlab="10個のサンプルの平均値", 
ylab="相対度数", main="平均値のばらつき")

cent

同時に描いた曲線が平均30,標準偏差ルート10(分散10)の正規分布だ.確かによく重なっている.

さて,次は分散比を調べよう.

まず,先ほどと同じように10個のサンプルを採取する.そして,もう一組10個のサンプルを採取する.二組のサンプル集団ができたら, 両方のサンプルそれぞれについて分散を計算する.そして一つ目のサンプルの分散を,二つ目のサンプルの分散で割る.

要するに,同じ集団から2組のサンプリングをしたときに,その「ばらつきの仕方」 の食い違いはどのようにばらつくのかを計算するわけだ.

これを利用するとA社の製品とB社の製品の間に「当たり外れの大きさ」の違いがあるのか,などを調べることができる(関連: 分散分析).

さて,10個のサンプリングを2組作り, それぞれ分散を計算して比率を計算してなんてことを何度も何度もやるのは気が遠くなるのでRに丸投げする.

sample.var <- numeric(0)
for(i in 1:2000){
 x <- sample(norm.tips, 10, replace=T)
 y <- sample(norm.tips, 10, replace=T)
 sample.var[i] <- var(x)/var(y)
 }

これで,"sample.var"という入れ物の中に2000個の分散比が入った.今回は早速ヒストグラムを描く.

par(mar=c(4,4,2,2))
truehist(sample.var,xlim=c(0,10),ylim=c(0,0.8),ylab="",xlab="")
par(new=T)
par(bty="n")
curve(df(x,9,9),xlim=c(0,10),ylim=c(0,0.8),ylab="相対度数", xlab="分散比")

varratio

このように,分散比は右側に長い尻尾を持つ分布をする.

ところで,このグラフに描いてある曲線は何なのか.

これは「F分布」と呼ばれる分布型の一つで,分散比のばらつきの仕方を表す.

先ほど,サンプリングの平均の分布の仕方(というか分散の大きさ)を決めたのはサンプリング数であった. ならばF分布もサンプリング数でその形が決まるのか.決まるのだ.

ただし,今回はサンプリングが2組ある.つまり,サンプリング数は2つある.なので, F分布の形は2つのサンプリング数によって決定される.より正確には,2つの「自由度」がその形を決定する.上の図に描いたF分布の曲線は 「第一自由度9,第二自由度9のF分布」と呼ばれる.「自由度」というのはこの場合「サンプル数-1」だが,必ずしもそうとは限らない. 詳しい説明はできないけど,たいていの場合「データ数-用いた平均値の数」が自由度となる. 上の例ではサンプルから分散を計算するときに平均値を一つ使っているため,「サンプル数-1」が自由度となる.らしい. 自由度に関する軽くて直感的に分かりやすい解説は「推計学のすすめ」を参考に,詳しい解説は…どこにあるんだろ?

ところでグラフの描画はpngとかepsとかのデバイスを直接開いて行ったほうが楽だということにいまさら気づいた.

ウィンドウを開いておいて微調整して,最終的に出来上がったものを保存する.というのはExcel的でまあ楽なんだけど, 時々プツっと落ちるのにいい加減嫌気が差していた.デバイスを開いておいてそこに描画していくようにすれば, スクリプトを流用することでデータだけ異なった体裁の全く同じグラフが何枚でも書ける.

また後でブログにやりかたをメモっておこう.今はデータがマジなデータしか無いからだめだ.

必要性に駆られるといろんなこと覚えられて楽しいな.疲れるけど.

posted by Rion778 at 00:13 | Comment(0) | TrackBack(0) | 勉強ノート | このブログの読者になる | 更新情報をチェックする
2008年01月25日

日本語,というか言葉が難しい

それ自身に関係する関係が他者によって措定されたのである場合には,その関係はもちろん第三者ではあるが,しかしこの関係, すなわち第三者は,やはりまたひとつの関係であって,その関係全体を措定したものに関係しているのである.
セーレン・キルケゴール 「死にいたる病」

まーちょっと落ち着きたまえと言いたくなるような文章だが, 各国語の訳も解説もあるのだから少なくとも誰かには伝わっているはずと考えて差し支えないだろう.僕にはイマイチ伝わってない. 前後の文脈を読むとぜんぜん分からないからちっとも分からないになり, 解説を読むことでかろうじて分かったよーな分かってないよーな感じを受ける.

思うに,材料が足りない.いや,誤解されるのを避けて材料をあえて使っていないのかもしれない. 「それ自身に関係する関係」を「自己」,「他者」を「神」と置けば少しは読みやすくなるだろうけど, それらの言葉に大なり小なり意図と異なりかねないイメージが付いていることは言うまでも無い.正確に伝えようとする故か, あるいは伝える人間に対する篩か.

どんなに奇抜なことを伝えようにも,手元に無い道具は仕えない.既存のものから作り上げる必要がある.

どんな建造物を建てるとしても,君はこれらの材料を,しかもこれらの材料だけを, やりくりして組み立てていかなければならない.
ウィトゲンシュタイン「論理哲学論考」

でも言葉というのはそのように出来ているものだ.つまり,言葉を使って伝えるならば常にそれは新しいものだろう. 常に既存の道具を用いて,新しいことを伝える.ただ,「既存」と「新しいもの」の関係は幾何学のように厳密ではなく,突然新しい (新しく見える)道具が出てきたり,平気で循環論法が使われたり,矛盾,誤用は時に「正解」となる.そこで重要なのは厳密さではなく, 相手の判断.言語を第一に形作るのは論理ではなく,他人だ.

その相手,その時間においてどのように理解されるかを考慮しておけば,言語を使うことはそれほど難しくは無い. 問題はこれらがいずれも予測の難しいものであるということだけだ.我々が言語を上手く使えないのは,単に状況を理解できないから. 我々が言語を理解できないのは,発信者の想定する状況を理解できないから.

問題は単純だが深刻だ.必要なのは譲歩か,更なる説明か.

posted by Rion778 at 01:38 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月22日

LaTeX関係の覚書(図表の代わりの空白,数字なしの章見出し)

図や表の代わりにとりあえず空白を入れておく.以下は図の場合.表の場合は\figureを\tableへ. キャプションは上につけたほうがいい.

\begin{figure}
  \begin{center}
          \vspace{2in}
          \caption{かくかくしかじかの図}
          \label{hoge}
  \end{center}
\end{figure}

章見出しとかの番号は消して,目次には章タイトルを出力する.以下は章見出しの場合. 一行目は普通の章見出しにアステリスクを付けたもの.この命令だけだとページには数字なしで出力されるが目次に出力されない. 二行目で目次へ出力する.

\chapter*{ほげほげ}
\addcontentsline{toc}{chapter}{ほげほげ}

人に見せただけでテキストは壊れる.完璧へはただ漸近のみが許される.だから現在地と完璧の間にはいつだって点を取ることができる.

無限に修正を続けるのであれば完璧と等号で繋ぐことは許されるだろうか?

posted by Rion778 at 20:27 | Comment(0) | TrackBack(0) | PC関連。HTMLとか,Linuxとか,Rとか | このブログの読者になる | 更新情報をチェックする
2008年01月19日

(笑)

必要とする予備知識の少ないテキストは優れたテキストであろうか.

ある言葉,概念が難しいのは,多くの場合それが複合的なものだからだ.A,Bという2つの概念があったとして, これらの共通部分が和集合より大きくなることはない.だから複合的な概念は必然的に適用範囲が狭くなり,それを使う機会も少なくなる. 使う機会の少ない概念は,難しいと感じてしまう.ただ,大抵の場合単に慣れてしまうことで問題は解決する. 複雑な概念を簡単な部品に分解せずに認識できるようになるために必要なのは,「慣れ」だ.

一旦それに慣れてしまえばわざわざ概念を部分へ分解することは全く無駄な作業となる.

今,一つ前の記事で私は「プリアンブル」という単語を使ったが, 数ヶ月前に初めてそれを見たときその意味をすぐに理解することはできなかった(preamble=序言という意味を知っていれば見当くらいついたかもしれない). プリアンブルを簡単に説明してしまえば,単に\begin{document}の前の部分というだけのことである.ならばなぜわざわざ 「プリアンブル」と呼ぶのかといえば,その意味が複合的なので「プリアンブル」という言葉に置き換えてしまったほうが都合がいいからだ. 実際,単に\begin{document}の前という説明のほか「文章が始まる前の部分」「文章の体裁を指定する部分」 「パッケージを読み込む部分」などなど「プリアンブルとはかくかくしかじかである」という説明は何通りもできる.つまり,「プリアンブル」 という単語を用いない場合,目的に応じて記述の仕方を変えなければならないということである.それでは都合が悪い. 意味が正確に伝わらない可能性も大きくなる.そこで,これらの意味を全て包含した「もの」として「プリアンブル」という単語を用いる. そうすれば,少なくともこの単語に慣れたひとにとって誤解の恐れはなくなる.

対象を「…はかくかくしかじかである」という形へ分解することは,同じように「…はかくかくしかじかである」 という分解ができる程度にその対象に慣れたひとにとっては単に読む労力と誤読の可能性を増やすだけである.

それでは,私は私のテキストを読むひとの予備知識を知ることができるかといえば,できない.できないが,予想はできる.特に, 読む人間の種類が限られているような場合,それはぼんやりと形を帯びてくる.境界は分からないが, だからといって領域が分からないということにはならない.視界の限界を見ることができなくても,視野の広さを知ることはできる.しかし, ぼんやりとでもそれを知るためにはその中に自分がいなければいけない.そしてその領域に「慣れて」いなければいけない.

つまり,勉強するだけではだめで,誰かに伝えなければいけない.伝わったのか,伝わらなかったのかを確認しなければいけない.

日常語となると事態はさらに深刻になる.

私は「(笑)」と「(笑」と「(藁)」と「(w」と「w」と「wwwww」の間になんとなく意味の違いを感じている.そして, それを了解しながらこれらを使い分ける人間がいることを認識している.一方,そんなことは気にも留めていない人間がいることも認識している. これらの表現の違いは,場合によっては単にニュアンスの違いのみならず,使用者の属す文化的背景すら浮き彫りにする. それでいながらそれらの人々を見分ける術は知らない.「(笑)」を基本とする(と私は認識している)これらの表現のうち, 一体どれを文章の末尾に付けるべきかという問題は,私にとって極めて深刻な問題なのだ.

や,まあ大抵の場合使ってみれば片付く問題なのだけれど.

posted by Rion778 at 23:45 | Comment(0) | TrackBack(0) | 考え事 | このブログの読者になる | 更新情報をチェックする

tabular環境の中で行を縦に結合[LaTeX]

multirowパッケージを使う.プリアンブルにて

\usepackage{multirow}

を忘れずに.

使い方は,結合する一番上の欄に「\multirow{結合する行数}{幅(指定せずに*入れておけばいいかも)}{内容}」. 結合されるほかの欄は空欄にしておく.

以下ソースのサンプルを.

\documentclass{jsarticle}
\usepackage{booktabs}
\usepackage{multirow}
\begin{document}
 \begin{tabular}{ccc}
  \toprule
          \multirow{2}{*}{縦結合}& 1 & 2 \\ \cline{2-3}
                                  & 3 & 4 \\
  \bottomrule
 \end{tabular}
\end{document}

booktabsパッケージは罫線で\toprule,\bottomruleを使ったので使用.好みの問題.\hline, \clineのみを用いるなら不要.

なお,罫線は結合している場所だろうがお構いナシに引かれる.なので, 結合した欄を除いた罫線が欲しいときは\cline{欄番号-欄番号}できちんと指定する.

出力結果はこんな感じ.

multirow

結合したい行数が偶数だったら役に立つけど,奇数行だったら要らないかな.真ん中に書いてあと空欄にすればいいだけだし.

posted by Rion778 at 20:47 | Comment(0) | TrackBack(0) | PC関連。HTMLとか,Linuxとか,Rとか | このブログの読者になる | 更新情報をチェックする
2008年01月17日

変換の度に劣化する。伝言ゲームのごとく。

やっぱりLaTeXは画像の取り扱いが難しい.

画像はPicasaで適当に編集した後,Inkscapeで文字やら何やら書き込んで, Inkscapeからビットマップでエクスポートする.そうするときれいなpngが出来上がる.なのでできればそいつをそのまま使いたい.

LaTeXのgraphicxパッケージはdvioutって指定してやればpngを扱えるはずなんだ.でもうまくいかない. dviまでは何のエラーも吐かずに出来上がるけど,プレビューしてみるとプラグインが見つからないとかで止まってしまう.何でかな.

jpgならいける.でも若干劣化してしまう.ペイントなんか使ってるからなのだろうか.

それでも,分割して編集してるほかのファイルからLaTeXを実行したりすると見れなくなる時がある.何でなんだ.

まあ隙間の空け方なら分かる.\vspace{2cm}とか図の代わりに書いておけばいい. \tableとか\figureとかで括っておけば参照番号もふれる.だから最悪糊とハサミがあればなんとかなるけど, やっぱここらで一度印刷してみないとダメだな.

posted by Rion778 at 02:11 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月14日

いつかやる≒いつまでもやらない

そしてもっとも深遠な問題が実はいささかも問題ではなかったということは、驚くべきことではない。

ほんの少し駒の進め方を変えただけで一節分のテキストが消えてしまった. 重要だと思っていた推定値を求める必要がなかった.ただそれだけ.よく確認すればすぐ分かっただろうに,無駄な時間を費やしてしまった.… 「無駄じゃない時間」があるのかと問われても答えられないけど.

このまま捨てるのももったいない気がするので,機を見てブログの肥やしにしよう.

内容は多分商の平均値と調和平均の関係の話.それ自体は何の役にも立たないかもしれないけど, 手に負えない推定値を求める練習くらいにはなるかもしれない.

まあ,いつかヒマなときに.…いつになるやら.実際僕はいつでもヒマなので,「やらない」 宣言に非常に近い気がしないでもないでも.

しかしこれで大幅に書くべき内容が減ったはずなのに手元のテキストはいまだ際限なく増え続け, さっぱり終わる気がしない.余分なものは限界まで削ってるつもりなんだけど.間に合うのかこれ?

posted by Rion778 at 21:39 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月11日

一日座ってただけなのにこの足の疲労は何だ

グラフが際限なく増えていく…

しまいにゃLaTeXさんに「Too many unprocessed floats.」とか怒られる始末.

ちなみに,同じ症状で困ってる人は図表の合間に適宜 \clearpage というコマンドを入れるといい.\clearpage と書いた部分までで一旦ページを作り上げてしまうので,「フロート多すぎ!」とか怒られなくて済む.それと, 上手く使えば配置の制御にも役立つ.

LaTeXはテキスト関連はほぼ完璧に仕上がるけどちょっと画像や表の配置が上手くない気がするな. まあこっちで指示してやれば済む話なんだけど.

しかしデータってのはどの程度まで提示すべきなんだろうか.言及しないもの,する意味のないものは当然載せてないが,それでも多い. 先例はもっと多いが,そんな理由で納得したくない.大体,そいつをエレガントだとは感じなかった.尊敬できないものに倣うことができようか.

それと,統計処理なんてする意味あるんだろうか.いや,統計処理と括ってしまっては言いすぎか.検定だ.棄却検定が分からない.

平均値や分散といった統計量を計算することには当然意味があるが,棄却検定はまだ意味が把握しきれない. 見れば分かるような差があるデータは大抵統計的にも差が出てくるし, 見て分からないような差のデータなら統計的に差があろうが役に立たないことが多い.「見れば分かる」データは大抵の場合最終的に「見て」 判断するものだろう.明らかに差のあるグラフに申し訳程度に「異なるアルファベット間には…」と書いてあるのを見るにつけ複雑な気分になる. グラフで示すなら標準誤差を付けておけば十分だろう.あるいは箱ヒゲで示せ.どうせ統計処理をするなら全てのP値を書けよと. P値が5%以下だと有意だとか,フィッシャー先生が目安に考えた値にいつまで縛られているんだ. 今はP値が簡単に計算できるというのにアスタリスクやらアルファベットだけで済ますなんて時代錯誤も甚だしい. 中途半端に検定するくらいなら基本統計量をしっかり載せておいて欲しい.

しかし視覚的に理解しやすいデータが得られる方法は有用な気もする.Tukeyの方法による平均値の差の同時信頼区間とか. ノンパラじゃないのと例数が増えると見易くなくなるのが難点だが,確率と程度が同時に一目の元に理解できるのはいい. それでもやっぱ生データを比較して得られるのと情報量はそれほど変わらない気もする.

分布型が想定できないとか,直感的な理解が困難であるとかいった場合以外は統計なんて役立たないんじゃないかと思う今日この頃. 直感的に理解しやすいデータは直感的に理解すればいいと思う.でもそんなケースはそれほど多くないか. 例えば水準が4つあるだけでも散布図4次元になっちゃうし.やっぱ要るか.統計は.でも棄却検定の意味はやっぱわからん.

posted by Rion778 at 00:00 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月10日

データのまとめ方が分からない

どうも比率という形になっているデータは扱いにくくて困る.

比率が,二項分布の形で表されるような,すなわち,1と0の形に変換されたデータの平均値であるならばまだいい.平均値を望むならば, 重み付け平均値を計算すればよい.1に近い,あるいは0に近い値を含むデータについて統計処理を望むならば, 逆正弦変換を施せばより適切な結果が得られる.Yes,Noで形成される「比率」の取り扱いに関する情報ならあふれている.

が,同じ次元を持つ連続量から形成される比率の扱いはどうすればいいのだろう.例えば,長さ/長さ,重さ/重さ,体積/体積などなど. 分母が同じなら良い.それは通常の測定値と同様の取り扱いをして何ら問題を生じない.しかし,分母が異なった場合,分母の「重み」 を考慮するか否かにより「平均値」は異なった値を示す.

例えば,一つ目のデータが「0.1/0.1」で2つ目のデータが「0/1」だったら,その「比率の平均」はいくつになるのだろう. 50%か?それとも約9%か?二項分布するデータと同様の扱いをすべきだと言うならば,後者が正しいことになる.しかし,この「データ」 はあくまで測定値だという立場をとるならば,前者が正しいことになりはしないだろうか.

もっとも,統計処理のみが目的ならば順位を基準とした比較をすれば問題はないだろう.

しかし,時にデータというのは視覚的にまとめられた状態での提示を求められる.そこで要求されるのは平均であり,分散であり, データ数である.一体どのようにまとめれば良いのだろうか.

はてブの人気エントリに今年も「東大で学んだ卒業論文の書き方」 が上がってきている.もうそんな時期なんだ.我々はそんな時期の真っ只中にいるのだ.

posted by Rion778 at 00:19 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月07日

雑記

エコの論文作法はたしかに役に立つ本だと思う.

僕が文系だったなら.

LaTeXで組むとどんな文章もそれっぽくなってしまうのは良いな. 内容が薄くても見た目はちゃんとしたものが出てくるからテンションが上がる. 参考文献のところだけ見てたらとても僕が作ったとはとても思えない(まあBibTeX使ってんだけど). プログラムのようにコメントを入れられるのも便利だ.文献番号を振れ,表を入れろ,図を入れろ,明日の自分へどんな指示も送れる. 目次だ文献番号だ式番号だといった細かい数字は勝手に処理してくれるから,僕はテキストだけに集中できる. 大学でLaTeX教えないってのはちょっとしたイジメだとすら思える.今から思うと, ワープロソフトで論文やレポートをきれいに仕上げるのには無駄な労力がかかりすぎる.

ところでヒマだったので「空気を読む」ことについてまだ考えていた.「常識的に考えて」の判断基準は何なのだろう. 実体はどこにあって,誰がくれるものなのだろうと.

実際,「一般的なもの」というのはその辺に転がっているからサンプルには事欠かない.けれど, そいつらは自分の内側に取り込まないと観察できない.そうなるともうそれは単に自分を観察しているにすぎないのではないか.それを「常識」 と呼んでいいのだろうか.結局のところ,「空気読む」ことや「普通にする」ことのためには「一般的なもの」をただ「集める」 というのが唯一の手段ではないか.いくら言語化しようとしたって無駄だろう.「普通って何なのか」と聞いてみたって,「普通は普通だよ」 としか「普通は」答えられない.トートロジーになってしまう.それは結局「論理」や「倫理」 がどうしてこういう構造をしているのかという話であって,そんな簡単に言語化できるものではないのかもしれない.「人は沈黙せねばならない」 ってこと?なんかやっぱりよくわからんな.

まあそんな話よりも手元の仕事は簡単なハズだからとっとと進めよう.

posted by Rion778 at 00:09 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月04日

completeという意味で

風邪で味蕾が…

これはちっとも楽しくない.

というか体重が思ったより減っている.風邪ひいてたけどヒマだったので料理はちゃんとしてたというのに. 10の位の値が動いたら死ぬような気がして仕方がない.

あーしかし書かなきゃならないテキストが書き進められない.

今まで自分は何をしていたんだろう?この調子で要らないものをどんどん切り落としていったら, 最後には何もなくなってしまうんじゃないだろうか.水を煮込んでも鍋には何も残らないわけだし.

もしも自分の仕事を誰かに伝えることが無かったら,その仕事はずっと「完璧」なままなのにな.ああ, それでも一応まとめる必要はあるか.でも,誰かに伝えたいっていう強い動機がないと苦労してまでテキストに起こそうなんて思わない. そういうモチベーションはどこに置いてあるんだろう.

んなこと言ってても仕方ないからデータの整理からはじめるか.

しかし,ヒマにまかせて手を出した余計な部分に自作の推定量とか入ってるんだけど扱いはどうすりゃいいんだろう. 推定値として適切なことはシミュレーションでしか確認していない.といっても誤差の大きさと標本分散, サンプルサイズの関係まではなんとなく分かっている.なんとなく,というのが曲者だが.

モンテカルロの結果を貼り付けておけば良いんだろうか.こういう場合はソースも貼らないといけないんだろうか. hogehogeしてて見苦しいから人に見せたくないな.

最尤推定量であることを示せばいいんだろうけど,どうすりゃ「示した」ことになるのか.尤度関数を微分?なにそれうまいの? って感じだしどうすりゃいいものか.

というかそもそも「おまけ」なのにこのままだと本筋より長くなるわ.

筋書きすらまともにまとめられない.憂鬱でしかたないな.

posted by Rion778 at 23:16 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする
2008年01月02日

J.C. Maxwell 「気体の動力学的理論の例示」を読んでみる

風邪で三半規管が…

…フワフワしてちょっと楽しいかも。

気体分子運動論をカジろうとMaxwell先生らの論文集を調子こいて買ってはみたものの、「気体の動力学的理論の例示 ( Illustrations of the Dynamical Theory of Gases )」という論文の 3ページ目でストップしてからすでに1ヶ月が経過している。別に放置しているわけじゃなくて、結構頻繁に考えてこの様。 もうちょっとマトモな頭が欲しいよ。

とりあえず今日は何もする気が起こらなかったのでこの部分を一日中考え続けた。 その結果どうにか読めたような気がしないでもないでもないでもないのでメモっておこう。

ここでの大筋の議論は「完全弾性球の運動と衝突について」。

気体分子を完全弾性球に見立てることで、その運動を記述してやろうというわけです。読んだのは、反跳の方向の確率の部分まで。

要するに、2個の球がぶつかったとき、跳ね返る方向とその確率の関係はどうなってんの?という話。

で、本文。

第一の命題は、「自身の質量に逆比例する速度で逆方向に運動している2つの球が互に衝突するとき, 衝突後のそれ等の運動を決定すること」とある。

まずは次のようなお絵かきをする。

20080102fig1

AP、BQを衝突前の、Pa、Qbを衝突後の運動方向と大きさとすると、中心線Nに直角な方向の運動量は変わらず、 垂直な方向の運動量は反転することが分かる。「質量に逆比例する速度」というのがポイントだな。どの粒子も運動量が等しいと。

まあここまではいい。見れば分かる。

問題は次の命題。

「衝突後の速度の方向が与えられた範囲間にある確率を見出すこと」

まず、衝突の起こる条件を考える。

衝突がおこるためには一方の球の進路が他方の球の進路に重ならなければならない。QP間の距離をSとすると、 「片方の球の進路を中心とする半径Sの円の中にもう片方の球の進路が納まっていなければならない」、というふうに言い換えることができる。

そして確率。

まず、片方の球の進路を中心とする半径Sの円の中のある範囲内にもう片方の球の進路が存在している確率を求める。「ある範囲」 としてはドーナツのような輪っかを考える。輪の外周までの距離を(r+dr)、内周までの距離をrとし、 このドーナツの面積の中に進路が存在する確率を算出しよう。

要するに(ドーナツの面積)/(円の全体の面積)を算出すればよいので、

20080102eq1

が求める確率となる。πを約分して消した後、分子を展開して少し整理すると

20080102eq2

ここで、drは(多分)微小量を表している。なので、2次の微小量であるdr2は無視して (ほんとに無視していいのか?と思うところだけど、式中の他の値に対して微小な量であるdrの変化に対してdr2は 「微小」なので、式中の値と比較すれば0としてもほとんど問題にならない、ということらしい)

20080102eq3

が目的の確率となる。

ところで、最初の図において∠APaをφとおくと、∠APNは1/2φであり、rはS sin 1/2φと表すことができる。そして、 drはS sin1/2(φ+dφ)からS sin1/2φを引いたものとして表すことができる。 テキストではちょっと見づらいのでTeXで書き直すと

20080102eq4

これを先ほどの式へ代入してみる。分母のS2がrとdrのSによって約分されるのはすぐ分かるので、次の形に。

20080102eq5

括弧の中、sin1/2(φ+dφ)を加法定理により展開

20080102eq6

括弧の中、sin1/2φの項を、-2sin1/2φで括る(多分ポイント)

20080102eq7

ここで今できた括弧の中に注目する。ここには、三角関数に関する公式の一つである「半角の公式」

20080102eq8

が使える。よって、

20080102ep9

ここで今sin2(1/4dφ)の項というのはdφに関する2次の項である(角度0付近で角度φの変化とsin φの変化は等しい)ため、先ほどと同じように無視し、

20080102eq10

が得られる。ここで、sin 1/2φ cos 1/2φは倍角公式より1/2 sinφであること、sin 1/2dφは1/2dφでほとんど近似できることを考えれば、

20080102eq11

…長かった。ていうか本当に合ってるだろうなこれ。論文にはこの式が載ってるけど、 sinがどこまで掛かってるかイマイチ分かりにくい(ちなみにさっきの20080102eq3から2行おいてこいつが出てくる)。

それじゃあ次にこの確率に対応する反跳の方向を考えよう。

反跳の方向は3次元空間上のあらゆる方向が考えられる。そういうときにある「範囲」を指定するには、 半径1の球の面積に対応させるのが便利だ。立体角を指定するときのステラジアンという単位と同じように。

先ほど求めた「確率」は反跳の方向がφからφ+dφの間にある「確率」である。なので次に求めるのは、「半径1の球において、 角度φからφ+dφの間にある球帯の面積」である。最終的には、「面積あたりの確率」を求めようとしている、ハズ、多分。

「球帯の面積」というのは、平行な2平面で球体を切り取ったとき、2平面にはさまれている部分の球の面積のこと。そして、 この面積は球の半径をr、球帯の幅(平面に垂直な方向)をhとすると、2πrhという非常に簡単な形で表される。 どの部分を切り取ったかは問題にならず、どんな幅で切り取られたかのみが問題となるわけだ(参考:球帯と球冠私的数学塾))。

で、半径を1とするなら、球帯の面積は2πhとなるので高さのみを求めればいいことになる。ここで「高さ」 というのは先ほど確率を求めるときに使ったrと直角な方向である。つまり、「高さ」はcosを使って

20080102eq12

と表される。加法定理により第2項を展開すると

20080102eq13

ここで先ほどの計算を思い出すと、第1項、第2項を上手くまとめて適当な変形をほどこせば無視できる項となることが分かる。また、 sin dφはdφで近似できるとして書き直すと、

20080102eq14

これが球帯の「高さ」だ。よって球帯の面積は

20080102eq15

さて、単位面積あたりの確率を出そう。

20080102eq11

20080102eq15

で割ればいい。すると、

20080102eq16

ただし、これまでの計算は「非常に小さい角度の範囲」を対象としてきた。なので、ωを「半径1の球面上の任意の小さな面積」とし、

20080102eq17

と書き直しておこう。これは「半径1の球面上の任意の小さな面積」を反跳の方向が通る確率である。この値は角度φに依存していない。 すなわち、「全ての反跳の方向が等確率で起こりうる」。

長い長いながーい!

結論は簡単なのに。

しかしこの部分、元の論文(日本語訳)のテキストにして1ページ少々(A5)にすぎない。しかもこれ、 1859年の9月21にAberdeenの大英学術協会の会合で「読まれた」らしい。「Philosophical Magazine for January and july, (1860)」と表紙にあるから、 論文として掲載される前の話だろう。次元が3つか4つは違う感じだ。「読まれた」としてどんな反応があったんだろう。 パパっと検算してサッと流す感じだろうか。やっぱ本職の物理屋さんはすごいんだろうな。

ところで論文中には「Clausius教授」というような名前がたびたび出てくる。熱力学の基礎が築かれていった時代の話だ。 当時はまだ原子の実在すら確認されていない。なので、Maxwell先生も論文中で

粒子は堅い,球形の弾性体であると言う代りに, もし好むなら粒子はその作用がある小さな距離だけ離れたところ以外では感じないでその距離に来た時に急に非常に強い反発力として現れる様な力の中心であるという事もできるだろう.

と微妙な言い回しを使っている。結局は分子を「球」 に見立てて議論を展開しているというのに妙な一文だが当時は必要だったんだろう。

Maxwell先生以降に気体分子運動論を研究することになるBoltzmann先生は原子論を主張し、 反原子論の立場をとる学者ら(皮肉なことに「モル」という言葉を初めに使ったとされるOstwald先生もその一人だった) と激しい論争を繰り広げた。結局原子の存在が証明されたのは1909年のこと。… その3年前にBoltzmann先生は自ら命を絶っていたのだけれど。

…まーとりあえずここまで読んだ。でも合ってんのか本当に?

ああ、しかし2ページ後をめくってみると広義積分が…

一体いつ読みきれるのか…

参考

posted by Rion778 at 22:30 | Comment(0) | TrackBack(0) | 勉強ノート | このブログの読者になる | 更新情報をチェックする
2008年01月01日

あけおめ

雪だー!

風邪ひいてダルイから外出はやめようかと思ってたけど、楽しそうだから少し散歩してこようかな。

年の初めだからといってこれといってすることもないけれど、とりあえず目標は立てておこう。

今年の目標:以下の三つを読めるようにする

  1. 英語
  2. C
  3. 譜面

何?空気はどうしたって?

ところで大したことじゃないけど、「自分は空気が読めるかどうか」を自分自身で確認する方法は存在しないよね。

社会学者の宮台真司という人が名づけた「吉幾三問題」というものがある。

 吉幾三は青森から東京に出てきたときに、喫茶店で働いていたわけだ。ある時マスターがいなくなっちゃったんで、 一日中自分が店を任された。客がやってきて、メニューを見て「ウインナーコーヒーはありませんか」と言った。吉幾三ははたと考えて 「聞いたことがない、しかしウインナーコーヒーと言うくらいだから、ウインナーとコーヒーのことだろう」と思って、コーヒーと、 フライパンであぶったウインナーを出したら、客は何の疑いもなく食べて帰っていった。
宮台、東浩紀の〒(郵便) 本を語る その2

指摘してくれる誰かがいなければ自身の間違いを認識することは無い。「自分が正しいと思っていること」は 「自分が」正しいと思っているに過ぎない。いつ、どこで「間違いだ」という指摘をされるのか、間違いを悟るのか分からない。

誰かとすりあわせをすることで「正しさ」は共有されるんだろうか。でも、そこでの事実は「自分が」 誰かと正しさのすりあわせをした、ということだけであって、実際に「外」とシェアをしているという気はしないな。 最終的な判断が自分の内側にしか存在できない以上、ある程度の不確実性を取り払うことはできない。ある色に「赤」 という名前をつけることはできるだろう。でも、本当に同じものを見ていると誰が教えてくれる?

でも実際にはそんなことは少しも問題にならず、コミュニケーションは成立してしまう。 第三者から見れば明らかに奇妙な形であったとしても。

コミュニケーションはどうしたら成立するんだろう。「空気が読めた」ことは誰が教えてくれるんだろう。

posted by Rion778 at 15:23 | Comment(0) | TrackBack(0) | diary | このブログの読者になる | 更新情報をチェックする

広告


この広告は60日以上更新がないブログに表示がされております。

以下のいずれかの方法で非表示にすることが可能です。

・記事の投稿、編集をおこなう
・マイブログの【設定】 > 【広告設定】 より、「60日間更新が無い場合」 の 「広告を表示しない」にチェックを入れて保存する。


×

この広告は1年以上新しい記事の投稿がないブログに表示されております。