2017年6月統計検定2級 反省会会場
日曜に受けた統計検定2級の公式解答がでてきました。
http://www.toukei-kentei.jp/wp-content/uploads/ans2017j_grade2.pdf
結論から言うと、24/35 、率としては68％となりました。うーん微妙。よくよく見てみると基本的なところが間違っていて全然理解できないことがよくわかりました。。というわけで反省会！赤のところが間違っているところで、コメントをつけたしておきました。ただ見直してもよくわかってないようで、いろいろまずい部分が浮き彫りになっています。
[1] ④
Ⅰ ○ TV保有率はずば抜けていてPC保有率と被らない故
Ⅱ ○ PC保有率の第3四分位数12都道府県はDVD/BD保有率と被らない
Ⅲ × SP保有率とMP保有率はオーバーラップするため高いかどうかの判断がつかない
[2] ③ MPの最高値は75～80でb, SPの最高値は80~85でa,
[3] ③ 茨城の中央値は28841で足していくとD内
[4] ② 同じく長野の第1、第3四分位数の値まで足していくと、BとD
[5] ⑤ 北海道のA～Fまで足しても20％に満たない（エ）、秋田のA～Bを足して10％に満たない（ウ）
[6] ③ 1668/995 = 1.68, 68%増
[7] ②
Ⅰ × 季節性の傾向ならば前年同月比でもよくね？
Ⅱ ○ 図より拡大している
Ⅲ × 伸び率は訪問数の差分に影響しない
[8] ③ 自信なし！自己相関係数とコレログラムをお勉強
=> やはりまちがっていた！でも正解の絞り方がわからない。。。①なのかー①なんだなー。
[9] ④ 目分量、強い負の相関だけど、-0.9まではないか
[10] ② 目分量、ちょうど間とおてる
<span style="color: #ff0000">[11] ③ F値が、推定値/標準誤差で -6.27, 自由度はサンプル数25-1</span>
=> いやー、どうどうと間違った解答をｗ 自由度は 25-2 みたいです。。。おいうことでたぶん⑤
[12] ④
① × 制度は悪くなる
② × 迷った。各層内の散らばりを小さくする必要があるのか？
③ × 渋谷でやったら渋谷ピーポーの傾向が強くなる。無作為にはならない
④ ○ 非がない
⑤ × 渋谷ピーポーのお友達紹介したら渋谷ピーポーの傾向が強くなる。
[13] ⑤ 自信なし。④だったか。
=> ③でした。ｗで恥の上塗り。
b = (X-Y)/2 - (ε1 - ε2)/2 だから、σ**2/2, なのかなぁ。。
[14] ④ 超パニックになっていたが、0.1 * 0.4 だわね
[15] ③ 規格外率 0.1*0.4 + 0.05*0.3 + 0.02*0.3
[16] ② 条件付き確率 (0.1*0.4) / [15]
[17] ③ 5C4 (0.7)**4 (0.3)**1
=> ああ、間違えてもうた。そうか、4戦目までに名人が3勝1敗の場合、4C3 (0.7)**3 (0.3)**1 で5戦目に名人が一勝するために、0.7をさらにかけるのね。すると、0.28812で①でした！
[18] ② 多分間違えた⑤が正解か。6戦目までで引き分けている確率 6C3(0.7)**3 (0.3)**3 = 0.185<
=> やはり⑤でした。稼げた問題だけに超もったいない2問です。
[19] ④ 自信ない。複数の標本が合わさるやつ苦手ね。
[20] ④
Ⅰ ○ 足しても引いても平均0は0
Ⅱ ○ [19]が正しいこと前提だと真だが、果たして
Ⅲ × わからない。。。
=> Ⅲは○だったようです！これはどういうことだ。
[21] ③ 自信ない。正規分布で、0.05の上側は1.64となるので、カイ２乗だから、1.64の２乗つーことで選んでみました。
=> 完全に導き方ミスりましたね。カイ２乗分布表から0.5 自由値１の値をみつければ。3.84ですので④
[22] ⓶ [21]から、Wn/2 = 1.355, カイ２乗分布でそれを超えるのが、自由度７のとき、じゃあ nは１足して８！ それらしいあてずっぽだし、[21]の解答に依存する
=> ラッキーパンチ！一定確率であてずっぽが当たるのは皆平等に与えられた権利です。これはやっぱり分布表からでしょうね。2nが0.05を超える直前の値は7ですね。するとnは8になるようです。
[23] できなかった。とりあえず③と埋めてみた
=> やっぱりミスってた。んでまだ答えがわからんぞ。がんばろ。
[24] ⓸ 標準誤差 = 標準偏差/root(n)、正直パニクった
[25] ⑤ 計算して出てきたのがこれだけど、ミスってないか不安。1.64で算出
[26] ①
Ⅰ ○ 不偏・・・だよな？
Ⅱ × 許容範囲を広めるので、広くなる
Ⅲ × 1/root(3)倍になるはず。。。
[27] ② 自由度はそれぞれ１引く
[28] ① 27*49/100
[29] ② (2-1)*(2-1)
[30] ④ あれ、足せばいいんですよね・・・？
=> あちゃー、不偏分散は(n-1)で割らなきゃダメです！割ったら ①ですね。
[31] ⑤ 地域は４つなので、自由度は３、F値は平均平方の比率
[32] ② Prが0.0405で５％以下、なので棄却できない
=> これはやばい。検定の基本的な部分を間違えて覚えちまっている。5％以下なので、帰無仮説を棄却！しなきゃだめ。①
[33] ② 問１５は理解するのに時間がかかりパニクった。時間なくなるし。31.9 + 0.3 * 28 - 4.43
[34] ① Prが全部5%未満なので採用
=> 有意に「正」となるパラメータね。。。ひっかかった。β2は負のパラメータだわさ。
[35] ⑤
① × 博士号と理論研究の関係性はわからない
② × PhDAgeが少ないが最盛期も少なくなる
③ × 理論研究者のほうが４年半ほど早く最盛期に到達する
④ × t値の判定から、４年半の差は統計的に有意である
⑤ ○ だから、理論研究やれば４年半早く最盛期迎えられるの！
なんどもいうようですが、やばい。合格できるのかはなぞですが、基本的な部分がことごとく間違った理解をしていることが明らかになりました。[8][13][20][23]についてはきちんと解答が導けないありさま。受かるかまだわかりませんが、まだまだ精進が必要！
統計検定試験その後の散歩
前回のエントリの通り、本日統計検定２級をうけてきました。実は午後に準１級も申し込んでいたのですが、ほとんど勉強できなかったうえ、２級でぎりぎりの状態だったのでやめることにしました。
会場は文京区の中央大学キャンパス。そこの近くに文京シビックセンターによってから帰ることにしました。文京区の区役所なんですが実はここの上が無料の展望台なんです。東京２３区のど真ん中にあるため、２３区のランドマークがたくさん見渡せるのでおすすめ。
統計検定会場となった中央大学のキャンパスです。いい眺めでしょ。
シビックセンターをあとにしたあと、ふらふらと歩いていたら神社が。そこには金毘羅宮とかかれていました。そういえばブラタモリでやっていたことを思い出してお参りに。
そこから変なスイッチが入ってしまいました。よっしゃこのまま歩いて神田明神に行ってみようじゃないの！って感じで、水道橋駅から神田川をのぼって御茶ノ水、神田明神に到着。
なんか巫女さんとか宮司さんとかがいそいそと歩いていて、行事をやっている様子。境内の真ん中に藁のワッカができていて説明書きを読むとその中を通りながら８の字歩くとご利益があるようだけど、さすがに一人では恥ずかしいので却下。お参りだけして後にしました。
もうここまで来たら湯島もいっちゃうでしょ。一度はいったスイッチぶっ壊れちゃってるようです。
境内に入るとなにやら儀式をしている音が。そこでは神前の結婚式がとりおこなわれているようでした。音に合わせて２人の巫女さんが舞を踊っているようで、思わず前前前世を唄ってしまいそうなのをこらえながらお参りを済ませました。
いや、ここまでいろんな神様にお願いしたから統計検定も大丈夫でしょう！とルンルンで途中で見つけた定食屋にはいって、遅めの昼飯でラーメン・チャーハンセットを幸せにいただきました。きっと神様にも私の祈りが届いたのでしょうか、もしくは普段の行いかよかったのか、店をでたら激しい雨が降り始めてきてしまい、駅までダッシュしたのですがびしょ濡れになるはめに。。。きっと神様が僕の浅はかな信仰心とあるく不謹慎ぶりを見抜いて罰をあたえたもうたのでしょう。
統計検定の結果も思いやられます。
２０１７年６月統計検定２級 自己チェック
2017年統計検定2級を受けてきました。とりあえずわからなかったところも含めて自分の解答をさらします。決して模範解答ではありません、ていうか、間違いだらけだと思いますので、つっこみ大歓迎です。落ちてたら超恥ずかしいですが、公式回答がくるまで落ち着かないのでやってみました。
[1] ④
Ⅰ ○ TV保有率はずば抜けていてPC保有率と被らない故
Ⅱ ○ PC保有率の第3四分位数12都道府県はDVD/BD保有率と被らない
Ⅲ × SP保有率とMP保有率はオーバーラップするため高いかどうかの判断がつかない
[2] ③ MPの最高値は75～80でb, SPの最高値は80~85でa,
[3] ③ 茨城の中央値は28841で足していくとD内
[4] ② 同じく長野の第1、第3四分位数の値まで足していくと、BとD
[5] ⑤ 北海道のA～Fまで足しても20％に満たない（エ）、秋田のA～Bを足して10％に満たない（ウ）
[6] ③ 1668/995 = 1.68, 68%増
[7] ②
Ⅰ × 季節性の傾向ならば前年同月比でもよくね？
Ⅱ ○ 図より拡大している
Ⅲ × 伸び率は訪問数の差分に影響しない
[8] ③ 自信なし！自己相関係数とコレログラムをお勉強
[9] ④ 目分量、強い負の相関だけど、-0.9まではないか
[10] ② 目分量、ちょうど間とおてる
[11] ③ F値が、推定値/標準誤差で -6.27, 自由度はサンプル数25-1
[12] ④
① × 制度は悪くなる
② × 迷った。各層内の散らばりを小さくする必要があるのか？
③ × 渋谷でやったら渋谷ピーポーの傾向が強くなる。無作為にはならない
④ ○ 非がない
⑤ × 渋谷ピーポーのお友達紹介したら渋谷ピーポーの傾向が強くなる。
[13] ⑤ 自信なし。④だったか。
[14] ④ 超パニックになっていたが、0.1 * 0.4 だわね
[15] ③ 規格外率 0.1*0.4 + 0.05*0.3 + 0.02*0.3
[16] ② 条件付き確率 (0.1*0.4) / [15]
[17] ③ 5C4 (0.7)**4 (0.3)**1
[18] ② 多分間違えた⑤が正解か。6戦目までで引き分けている確率 6C3(0.7)**3 (0.3)**3 = 0.185
[19] ④ 自信ない。複数の標本が合わさるやつ苦手ね。
[20] ④
Ⅰ ○ 足しても引いても平均0は0
Ⅱ ○ [19]が正しいこと前提だと真だが、果たして
Ⅲ × わからない。。。
[21] ③ 自信ない。正規分布で、0.05の上側は1.64となるので、カイ２乗だから、1.64の２乗つーことで選んでみました。
[22] ⓶ [21]から、Wn/2 = 1.355, カイ２乗分布でそれを超えるのが、自由度７のとき、じゃあ nは１足して８！ それらしいあてずっぽだし、[21]の解答に依存する
[23] できなかった。とりあえず③と埋めてみた
[24] ⓸ 標準誤差 = 標準偏差/root(n)、正直パニクった
[25] ⑤ 計算して出てきたのがこれだけど、ミスってないか不安。1.64で算出
[26] ①
Ⅰ ○ 不偏・・・だよな？
Ⅱ × 許容範囲を広めるので、広くなる
Ⅲ × 1/root(3)倍になるはず。。。
[27] ② 自由度はそれぞれ１引く
[28] ① 27*49/100
[29] ② (2-1)*(2-1)
[30] ④ あれ、足せばいいんですよね・・・？
[31] ⑤ 地域は４つなので、自由度は３、F値は平均平方の比率
[32] ② Prが0.0405で５％以下、なので棄却できない
[33] ② 問１５は理解するのに時間がかかりパニクった。時間なくなるし。31.9 + 0.3 * 28 - 4.43
[34] ① Prが全部5%未満なので採用
[35] ⑤
① × 博士号と理論研究の関係性はわからない
② × PhDAgeが少ないが最盛期も少なくなる
③ × 理論研究者のほうが４年半ほど早く最盛期に到達する
④ × t値の判定から、４年半の差は統計的に有意である
⑤ ○ だから、理論研究やれば４年半早く最盛期迎えられるの！
あらためて見直してみると自信ないし基本的なところ抑えられてないし、めっちゃ恥ずかしい。自信ないとこ全部だめだったら落ちてるなぁ、こりゃ。答え出てきたら改めて採点してみます。
ランダム迷路をコーディング
モチベーション
昔よく遊んだゲームで「トルネコの大冒険」というのがありました。いわゆるローグ系の走りで毎回変わるダンジョンの配置に驚かされたものです。今回はそのダンジョンをランダムに作る・・・と思ったのですが、頭で考えた限りではちょっと設定がめんどくさそうかなぁと。そこで、トルネコのダンジョンのかなり下の階にでてくる、ランダム迷路をつくってみようかと思いました。いまの私の実力にちょうどいい。よし、それをPythonでやってみようぞ！
まずは結果から
こんなんを出力しました。にょきにょき伸びていく感じがなんとも気持ち悪くて俺好みですよ。
きちんと迷路になっています。定義としてはループをつくらない、マスを全部埋めていく、ということを念頭においてました。まあ、よくみると埋まっていないますもあるのですが、まあ、おおめに見てください。
全部のコード
import numpy as np import random import matplotlib.pyplot as plt import matplotlib.animation as animation import matplotlib.patches as patches %matplotlib inline from IPython.display import HTML class MakeMaize(): def __init__(self): self.size = [60,60] self.holes = [[[1,1,1,1] for i in range(self.size[1])] for j in range(self.size[0])] self.holes[0][0] = [0, 0, 1, 1] self.driller = [] self.died = [] self.dlog = [] #上下左右のマスに動けるかチェック。[下、左、上、右]でいけたら1 def check_available(self, plr): direction = [1,1,1,1] # 今来た道は戻れない direction[(plr["f"]-2)%4] = 0 # 壁にいるときは壁側に移動できない if plr["p"][0] == 0: direction[1] =0 if plr["p"][1] == 0: direction[0] =0 if plr["p"][0] == self.size[0] - 1: direction[3] = 0 if plr["p"][1] == self.size[1] - 1: direction[2] = 0 # それ以外で、すでに埋まっている場所を確認(holesに記録) x = plr["p"][0]; y = plr["p"][1]; for i, d in enumerate([[0,-1],[-1,0],[0,1],[1,0]]): dx = x + d[0]; dy = y + d[1]; if direction[i] ==1: #print "DX:DY:", dx, dy, self.holes[dx][dy] ix = 1 for j in self.holes[dx][dy]: ix *= j direction[i] = ix # 結果を返す return direction # ３から９までの値をランダムに返す。分岐、折れ曲がりのタイミングを設定 def ran_gene(self): return random.choice(range(3,10)) # ルーティン drillerが場所を定義 def step_forward(self, c): def direc(f, p, g): x, y = 0, 0 if f == 0: y = - 1 if f == 1: x = - 1 if f == 2: y = 1 if f == 3: x = 1 while len(self.dlog) != c+1: self.dlog.append([]) self.dlog[c].append([[p[0], p[0]+x], [p[1], p[1]+y]]) return [p[0]+x, p[1]+y] # Drillerの１ステップ for d, drl in enumerate(self.driller): # 現在の動ける場所を確認 av_dire = self.check_available(self.driller[d]) # 動ける場所がない場合は動ける場所にワープ if av_dire == [0, 0, 0, 0]: cands = [] # 過去ログ中からランダムに動ける場所に移動 pick = random.choice(self.dlog) random.shuffle(pick) for p in pick: cand = self.holes[p[0][0]][p[1][0]] if sum(cand) > 0: ordr = [0, 1, 2, 3] random.shuffle(ordr) for o in ordr: if cand[o] == 1: self.driller[d]["f"] = o self.driller[d]["p"] = [p[0][0], p[1][0]] break if self.driller[d]["p"] == [p[0][0], p[1][0]]: break # ターン、分岐のカウントを一つずつ減らす self.driller[d]["tc"] -= 1 self.driller[d]["sc"] -= 1 # 現在の進む向き f = self.driller[d]["f"] av_dire[(f-2)%4] = 0 # 現在の向きに進める場合、ターンするカウントが0出ない場合 # 現在の方向に一つ進む if av_dire[f] == 1 and self.driller[d]["tc"] != 0: self.driller[d]["p"] = direc(f, self.driller[d]["p"], d) av_dire[f] = 0 # 進める方向からランダムに一つ方向をえらぶ # 前に進めない場合の対処 way = [] for i, v in enumerate( av_dire): if v == 1: way.append(i) if len(way) > 0: af = random.choice(way) # ターンカウントが0のとき # 進める方向にターンする if self.driller[d]["tc"] <= 0: self.driller[d]["p"] = direc(af, self.driller[d]["p"], d) self.driller[d]["f"] = af # 分岐カウントが0のとき # 新しいdrillerを作成してリスト追加、進める方向に進ませる elif self.driller[d]["sc"] <= 0: #print "\tAF", af, av_dire, self.driller.append({ #"p": direc(f, self.driller[d]["p"]), "p": self.driller[d]["p"], "f": af, "sc": self.ran_gene(), "tc": self.ran_gene() }) av_dire[af] = 0 else: if d not in self.died: self.died.append(d) # 各カウントが0になったときにランダムに値を設定 if self.driller[d]["tc"] == 0: self.driller[d]["tc"] = self.ran_gene() if self.driller[d]["sc"] == 0: self.driller[d]["sc"] = self.ran_gene() self.holes[self.driller[d]["p"][0]][self.driller[d]["p"][1]] = list(av_dire) return 0 # ログの表示 def print_now(self, c): print c log = [] for i, drl in enumerate(self.driller): print " ", "PLY", i, "P:", drl["p"], "F:", drl["f"], "SC:", drl["sc"], "TC", drl["tc"] log.append(drl["p"]) self.logs.append(log) def start(self): # 初期のdrillerの配置 self.driller.append({"p": [0, 0], "f": 3, "tc": self.ran_gene(), "sc": self.ran_gene()}) #self.driller.append({"p": [self.size[0]-1, self.size[1]-1], "f": 1, "tc": self.ran_gene(), "sc": self.ran_gene()}) # 各drillerを動かす。 c = 0 while(True): #self.print_now(c) code = self.step_forward(c) if code < 0: break if len(self.died) == len(self.driller): break c += 1 if c == 500: break dlog = [] p_old = [] print len(mm.dlog)
これを使って可視化します。
mm = MakeMaize() mm.start() fig = plt.figure() ax = plt.gca(axisbg='black') ax.set_xlim(-1, mm.size[0]+1) ax.set_ylim(-1, mm.size[1]+1) fig.patch.set_facecolor('black') def init(): #line.set_data([0,0], [0,0]) #return line, return [] def animated(i): patches = [] for xy in mm.dlog[i]: x, y = xy patches.append(ax.add_patch( plt.Rectangle((x[0], y[0]), x[1]-x[0], y[1]-y[0], color="white", edgecolor="white", linewidth="2.5"))) return patches ani = animation.FuncAnimation(fig, animated, frames=len(mm.dlog), init_func=init, interval=100, blit=True) ani.save("makemaize.gif",writer='imagemagick') HTML(ani.to_html5_video())
嫁のいない土曜の昼は・・・
さぁ、土曜日だ。嫁も朝から外出中。こんな昼はそうだあれしかない。嫁がいない休みの昼にしか出ないことをするんだ。そう、それは・・・
そうめん祭りじゃ～
ふぅーーーーーー、いやっはっ！！！
そうだよ、そうめんをたくさん食うなんて暴挙は一人じゃなきゃできない暴挙、贅沢、（パチン、）すしざんまい！！！一束、二束、チっチっチっ、嫁がいない土曜の昼になんたる甘ったれたことを・・・。ふざけんなよ、三束に決まってるでしょーが！ふひゃーい、待ってろよ！
んっ、何、食生活の偏りが心配だと・・・？はっはっはっ、その点は抜かりはない。これをみよ！
とりごもくおむすびとポテトサラダ
もう一度いう、とりごもくおむすびに、ポテトサラダ。かつてこんな最強な組み合わせがあっただろうか、とりごもく、ポテトサラダ、もちろんメインはそうめんだっ！！！ナイスバランス！！SMAPでいえば中居くんのような絶っ妙なバランス感覚ぅ～。 向かうところ敵なし、鬼に金棒、さぁこれでみんなも恐れおののいたはずだっ！
さて、メインだよメイン。そうめんゆでちゃうよ。三束だよ、三束。はっはっはっ、
ここから先のお楽しみにはおめぇらにレポートするわけにはいかねぇ。なにがおきたかは内緒だ。これだけ伝えてやっただけでも感謝してほしい。それじゃぁまたPythonの世界で会おう！
ライフゲームをコーディング
モチベーション
あたらしい人工知能の教科書 プロダクト／サービス開発に必要な基礎知識
- 作者: 多田智史
- 出版社/メーカー: 翔泳社
- 発売日: 2017/01/10
- メディア: Kindle版
- この商品を含むブログを見る
この本読んでます。知っていることから知らないこと、新しい知識がいっぱいつまっていて読み応えは十分。 そのなかでライフゲームというゲーム？があることを初めて知りました。ルールは以下の通りです。
大きなグリッド内で一マスの生物が生きている、もしくは死んでいる状態で存在している。 隣り合う8つのマスの状態でそのマスの次の生存状態がきまる - 死んでいるマスの周りで3つの生きているマスがいると、次にそのますは生きている状態になる - 生きているマスの周りで生きているマスが1つ以下だと、次にそのマスは過疎で死ぬ - 生きているマスの周りで生きているますが4つ以上だと、次にそのますは過密で死ぬ
それ、Pythonでやってみましょう！
コーディング
import numpy as np import random import matplotlib.pyplot as plt import matplotlib.animation as animation import matplotlib.patches as patches %matplotlib inline from IPython.display import HTML class Seizon(): def __init__(self): self.h = 50 self.v = 50 self.space = [[0 for i in range(self.h)] for j in range(self.v)] self.cycles=15 self.start_val = 0 self.rectPtrn = [] def nextage(self, x, y): a_cell = 0 for i in [-1, 0, 1]: for j in [-1, 0, 1]: if x+i < 0 or x+i > self.h - 1 or \ y+j < 0 or y+j > self.v - 1 or (i==0 and j==0): pass elif self.space[x+i][y+j] == 1: a_cell += 1 if self.space[x][y]==0: if a_cell == 3: return 1 return 0 if self.space[x][y]==1: if a_cell == 2: return 1 return 0 def check_life(self): for i in range(self.v): for j in range(self.h): test = self.nextage(i,j) self.space[i][j] = test def random_set(self, d=2): c = 0 for i in range(self.v): for j in range(self.h): self.space[i][j] = random.randint(0,99) % 2 c += self.space[i][j] self.start_val = c def calc_alives(self): c = 0 for i in self.space: for j in i: c += j return c def pickRect(self): h = 0 v = 0 getRect = [] for i in self.space: v = 0 h += 1 for j in i: v += 1 if j == 1: getRect.append([v,h]) self.rectPtrn.append(getRect) def __main__(self): self.random_set() for t in range(self.cycles): val = self.calc_alives() print str(t) + " : trial ", print val * 1.0 /self.start_val, val self.pickRect() self.check_life() print len(self.rectPtrn) sz = Seizon() sz.__main__()
出力
0 : trial 1.0 1236 1 : trial 0.532362459547 658 2 : trial 0.366504854369 453 3 : trial 0.25 309 4 : trial 0.155339805825 192 5 : trial 0.0873786407767 108 6 : trial 0.0647249190939 80 7 : trial 0.0493527508091 61 8 : trial 0.0404530744337 50 9 : trial 0.0364077669903 45 10 : trial 0.037216828479 46 11 : trial 0.0380258899676 47 12 : trial 0.037216828479 46 13 : trial 0.0355987055016 44 14 : trial 0.0364077669903 45
ランダムで50%程度のマスカラ始めた場合、10回前後で大体3.5%で生き残り続けるようですね。この値は最初の配置で変わるようですが、3~4%ぐらいのようです。では実際どんな感じで生きているんでしょうか。アニメーションを使ってみてみましょう。
アニメーションをつかって可視化
fig = plt.figure() ax = plt.gca() ax.set_xlim(1,50) ax.set_ylim(1,50) #ax.set_aspect(1) def init(): return [] def animated(i): plt.cla() patches = [] print i, len(sz.rectPtrn[i]) for r in sz.rectPtrn[i]: #for r in [[1, 1], [2,2], [3,3]]: patches.append(ax.add_patch(plt.Rectangle((r[0],r[1]), 1, 1, fc="Blue"))) return patches ani = animation.FuncAnimation(fig, animated, frames=sz.cycles, init_func=init, interval=1000, blit=True) HTML(ani.to_html5_video())
こんなんが出ました。面白いですね。きまった型になってそのままキープしている様子がうかがえます。
天邪鬼的に修正天動説の動画GIFをコーディングしてみた
モチベーション
僕は生まれてこのかた、自分の感覚的に地動説が正しいと感じることがありません。この地球上にいる以上は天動説こそが真理に感じてしまいます。
ではなぜそれが正しくないということを知っているかというと、こちらのエントリにかかれているような科学の巨人達の肩の上にのって高等教育を受けさせてもらったことに他ならないのです。百聞は一見にしかずとは言いますが、目で見たことが正しいとも限らない、ということがあるのでなかなか難しいです。
一方で天動説を修正してみれば地球を中心とした正しい太陽系の図も作れるんですよね。ようはどこから観測するかの違いですもんね。地動説は神の視点でみたものですが、この修正天動説は地球上の観測者を固定して、より僕たちの実感に近いものになるはずです。
それ、Pythonでつくってみました
まずは地動説
ソースは後で掲載します。Jupyterコーディングしてアニメーションにしてみました。見やすさを考えて、太陽、水星、金星、地球、火星、月だけ。距離と速さは実際と同じ比率で設定していますが、月と地球の距離を実際比にするとあまりに地球に近くなるので実比よりだいぶ話しています。面倒なんで、月の公転速度も単純にして、太陽を回る地球の公転速度の1/12にしちゃいました。
まあ、見慣れていますよね。実にシンプル。念のため、真ん中が太陽で、内側から水、金、地、火の順に好転してます。地球の周りをまわっているのが月ですね。
続いて修正天動説
さて、地球を中心とするとどのようになるのでしょうか。。。
おお、なんかきしょくわるいっすね。真ん中が地球です。内側の円が月、外側が太陽になります。当然ほかの惑星は太陽の周りをまわっているんですね。
無理やりであることは承知の上ですが、こうした地球系のソーラーシステムというのも、閉じた宇宙空間では間違いではないと個人的に思うんですよね。もちろん、宇宙空間を説明するうえで、地動説がシンプルに働く、ということではあるのですが。
今回はJupyterの上で、Python2.7を使ってコーディングしました。まずは地動説のコード。超適当かつざっくり貼り付けます。
```python
import math
import matplotlib.pyplot as plt
import matplotlib.ticker
import matplotlib.animation as animation
%matplotlib inline
%matplotlib nbagg
class Planet():
def __init__(self, ax, sp, dm):
self.ax = ax
self.sp = sp
self.dm = dm
def diameter(self):
return self.dm
def axis(self, time=0):
if self.ax == 0:
return (0.0, 0.0)
else:
return Ear.position(time)
def position(self, time=0):
spd = 0.0
if self.sp != 0: spd = 1.0/self.sp
rd = time * 2*math.pi/10 * spd
if self.ax == 0:
return math.cos(rd) * self.dm, math.sin(rd) * self.dm
else:
x, y = Ear.position(time)
return x + math.cos(rd)* self.dm, y + math.sin(rd)*self.dm
Sun = Planet(0, 0.0, 0.0)
Mer = Planet(0, 0.24, 0.39)
Ven = Planet(0, 0.62, 0.72)
Ear = Planet(0, 1.0, 1.0)
Mar = Planet(0, 1.89, 1.5)
Moon = Planet(1, 0.083, 0.1)
plnt = [Sun, Mer, Ven, Ear, Mar, Moon]
fig = plt.figure()
ax = plt.axes(xlim=(-1.6, 1.6), ylim=(-1.6, 1.6))
ims =
for p in plnt:
if p.ax==0:
circle = plt.Circle( p.axis(0), p.diameter(), fill=False)
ax.add_artist(circle)
for i in range(300):
x, y = ,
for p in plnt:
px, py = p.position(i*0.1)
x.append(px)
y.append(py)
im = plt.plot(x, y, "bo")
ims.append(im)
ani = animation.ArtistAnimation(fig, ims, interval=100)
plt.show()
```
そして修正天動説。
```python
import math
import matplotlib.pyplot as plt
import matplotlib.ticker
import matplotlib.animation as animation
%matplotlib inline
%matplotlib nbagg
class Planet():
def __init__(self, ax, sp, dm):
self.ax = ax
self.sp = sp
self.dm = dm
def diameter(self):
return self.dm
def axis(self, time=0):
if self.ax == 0:
return (0.0, 0.0)
elif self.ax == 1:
return Sun.position(time)
else:
return Ear.position(time)
def position(self, time=0):
spd = 0.0
if self.sp != 0: spd = 1.0/self.sp
rd = time * 2*math.pi/10 * spd
if self.ax == 0:
return math.cos(rd) * self.dm, math.sin(rd) * self.dm
elif self.ax == 1:
x, y = Sun.position(time)
return x + math.cos(rd)* self.dm, y + math.sin(rd)*self.dm
else:
x, y = Ear.position(time)
return x + math.cos(rd)* self.dm, y + math.sin(rd)*self.dm
Sun = Planet(0, 1.0, 1.0)
Mer = Planet(1, 0.24, 0.39)
Ven = Planet(1, 0.62, 0.72)
Ear = Planet(0, 0.0, 0.0)
Mar = Planet(1, 1.89, 1.5)
Moon = Planet(0, 0.083, 0.1)
plnt = [Sun, Mer, Ven, Ear, Mar, Moon]
fig = plt.figure()
ax = plt.axes(xlim=(-2.0, 2.0), ylim=(-2.0, 2.0))
ims =
for p in plnt:
if p.ax==0:
circle = plt.Circle( p.axis(0), p.diameter(), fill=False)
ax.add_artist(circle)
for i in range(300):
x, y = ,
for p in plnt:
px, py = p.position(i*0.1)
x.append(px)
y.append(py)
im = plt.plot(x, y, "bo")
ims.append(im)
ani = animation.ArtistAnimation(fig, ims, interval=100)
plt.show()
```
その他メモ
動画GIFをつくるのに以下のサイトを参考にしました。
とても苦労したのは ffmpeg や ImageMagick のファイルパスを .matplotlibrc の設定ファイルに追加したのですが、クォーテーションでくくっているとダメなんですね。これを理解するのに6時間がかりでした。あー、無駄な時間を。。。
