情報地球科学演習 (2024年前学期) 第05回

前回までの補足

「gfortran や Linux が怒っているのは分かるのですが、何に怒っているのかが分かりません。」

前回 にも指摘した通り、まずは 冷静にエラーメッセージを読んでみる べし。 一見すると「冷たい」「不親切」のようだが、実はどう直せばいいか? のヒントをちゃんと教えてくれている。

コンパイル時のエラー

これまでの演習でよく見たコンパイル時のエラーメッセージと、その対処法は以下の通り。 gfortran が 何行目でエラーを出したか を確認して、その行に以下のような問題点がないかを調べてみよ。

Error: Symbol 'hogehoge' at (1) has no IMPLICIT type

hogehoge という変数が宣言されていない。 ファイルの先頭に戻って、変数 hogehoge を宣言せよ。 さもなければ、何か別の名前の変数を間違って hogehoge とタイプしていないか確認せよ。

Error: Symbol 'hogehoge' at (1) already has basic type of REAL

hogehoge という名前の変数が2回以上宣言されている。 不要な宣言を消去せよ。

Error: Unclassifiable statement at (1)

ふつうの整数とか実数とか、「配列ではない」ものとして宣言されている変数を、「配列」のように添え字をつけて使っていないか? その行の中に含まれている変数がどう宣言されているか を一通り確認してみよ。

Error: Unexpected STATEMENT FUNCTION statement at (1)

「配列」のように「番地つき」で参照しようとしている変数なのに、実際にはその変数が宣言されていないのではないか? 「配列」として使っている変数の名前をタイプミスしていないか 確認せよ。 ちなみにこのエラーメッセージそのものの意味は、(例えば sqrt(n) みたいな「カッコつき」だから) gfortran は「関数」だと判断したのに、実際には「関数」ではなかった、と言っている。しかし亀山がこのメッセージを見てきた例のほとんどは、「配列」として宣言した変数のタイプミスが原因。

Error: Incompatible ranks 0 and 1 in assignment at (1)

代入文の「=」の右と左で、配列とそうでない変数が本当に正しく区別されているか? 例えば、「配列ではない」ふつうの変数に、「配列」そのものを代入 しようとしていないか?

Error: Invalid character in name at (1)

何か不用な文字が入っている。 例えば「行末なのにコンマが入っている」とか。

Expected a right parenthesis in expression at (1)

最後にあるべき右カッコ「 ) 」がない。 エラーが出てきた式を見直し、左カッコ「 ( 」と右カッコの数や対応関係が間違っていないか確認せよ。 それでもまだエラーが出る場合には、1行の最大文字数 (132文字) を超えてしまっていないか確認せよ。

Error: Expecting END DO statement at (1)

end do の数が足りない。 do と end do の対応が正しくついているかを調べ、適切なところに end do を追加せよ。

Error: Expecting END PROGRAM statement at (1)

上の反対で、end do の数が多い。 do と end do の対応を見直し、不要な end do を消去せよ。

このうち、「end do の数」が多いとか少ないとか、あるいは「カッコの対応関係」が正しくないとかいう失敗は、前回 に紹介した emacs の例の「お役立ち」機能 を活用していれば自然と避けられるはず。 また、バグ取りの目的などでプログラムに手を入れる際には、いきなり疑わしい箇所を「ばっさり消す」のではなく、まずは 行の先頭に「!」をつけてコメント化 するのがよい。 その「ばっさり消し」てしまった部分に間違いがあるとは限らないのだから。

実行時のトラブル

これまでの演習でよく見た実行時のトラブルとその対処法のは以下の通り。 そのほとんどは、計算結果の出力に関するものだった。 書式仕様の詳しい説明は ここ を参照のこと。

Fortran runtime error: Expected INTEGER for item 15 in formatted transfer, got REAL

出力させようとするデータの数や変数型 (整数か実数か) が、書式仕様と一致していない。 1行にやたらと多くの値を書き出させようとしてはいないか?

計算結果に 「*」 がやたらと出てくる。

もちろん答が間違っているからなのだが.... 書き出す際の書式仕様を write(*,*) のように変更してもう一度コンパイル・実行し、どのくらい間違った答になっているかを調べてみよ。 なお、write(*,*) という出力書式の指定とは、「小数点以下の表示に何桁使うか?」などを gfortran におまかせ にすることの意味。

計算結果が 「NaN」 とか 「Inf」 になった。

多くの場合は、「負の数の平方根をとろうとした」 とか、「0 で割り算をしようとした」 というような、数学的に正しくない計算をさせようとしている ことが原因。 根号の中や割り算の分母に登場している変数の値が正しく計算されているかを確認してみよ。 ちなみに NaN とは Not A Number (非数)、Inf とは Infinite (無限大) の意味。

プログラムの実行が終わらない

どこかに 無限ループ を仕込んでしまったのではないか? ひとまず Ctrl-C (Ctrlキーを押しながら c キーも押す) してプログラムの実行を中止させ、プログラムの中身をもう一度チェックしてみよ。 例えば 汎用 do ループの中に「ループを抜ける」設定が抜けていたりはしないか?

あるいは、実に単純なことだが

$ cat matsuyama-kion.txt | ./a.out

として入力データをプログラムに与えないといけないところを

$ ./a.out

としてしまっただけというオチではないか。

実行するタイミングによって出てくる答が違う or 配列のある一部分だけ値がおかしい

値を正しく「初期化」しないままで変数を使っているときによく出てくる症状。 例えば「総和を入れる」用の変数の「中身をまず空っぽにする」操作を忘れていないか?

gfortran に見逃がされてしまうバグを見つける有効な方法は、プログラムの動作の途中経過を追跡 することであり、具体的には 「プログラムのこの時点ではこの変数にどんな値が入っているか?」 をつぶさに調べるのがよい。 「原始的な方法だなぁ」だと思うかも知れないが、プログラムの動作の不具合を引き起こすバグを見つけ出そうとするなら、これ以上に効果的な方法はない (残念だが本当)。 実際やってみると、プログラム本体を「気合いで」ひたすら読み返すのと比べて、こっちのほうが格段にバグを見つけやすくなることが分かるだろう。 亀山も本当に何度も体験しているのだが、「この俺様がそんな低レベルな間違いをするはずがない!!」とかいう思い込みが (やっぱり) あるせいで、プログラムをいくら「慎重に」見直したつもりでも、それだけではなかなか間違いを発見できないものなのです....

今日のテーマ: 最小二乗法

このテーマでの目標は、最小二乗法により「1月から12月までの各月の平均気温が過去134年間にどういう変化をしたか」をうまく「あてはめ」る直線の方程式を求め、さらにその「あてはめ」の決定係数の値も求めること。 ただし、今日の課題をスムーズに進めるには、前回 のゴールである「1月〜12月までの各月の気温の平均・分散・標準偏差を求める」プログラム (enshu03c.f90という名前だったか) が既に完成済であることが求められる。 このプログラムがまだできていない人は、今日のテーマに移る前に必ず仕上げておくこと。

最小二乗法による直線への「あてはめ」 (その1)

最初の目標として、「過去134年間の5月の平均気温の変化を直線で『あてはめ』る」ことを考えよう。 その際、グラフの横軸 (x にあたるもの) には西暦、縦軸 (y にあたるもの) はその月の平均気温をとることにする。

最小二乗法により2つのデータの間の関係を表わす直線の方程式 ( \( y=ax+b \) と書こう) を決めるには、その直線の 直線の傾き (\(a\)) と y 切片 (\(b\)) を求める必要がある。 それには講義で示した通り、既に求めた y の平均値の他に

• x の平均と分散
• x と y の積の平均

を新たに求めておく必要がある。 とはいえ、これらを求めるにあたっても、基本的な考え方はこれまでと全く同じである。

今日の演習でも、前回に作成したプログラム [みほん] の 適切な場所に、以下に示す「部品」をうまく追加する ことで目的のプログラムが完成するようになっている。 これには例えば

$ cp enshu03c.f90 enshu04a.f90
$ gedit enshu04a.f90 &

として「1月〜12月までの各月の気温の平均・分散・標準偏差を求める」プログラムを enshu04a.f90 というファイルにコピーし、このファイルを編集していくのが安全。 ただし、どこに追加するのが適切かは 各自で考える こと。 なお以下では例の如く、薄い文字で書いてある行は、新たに書き足す必要のない 行であることを示す。 ということは、これまでに作ったプログラムの中にこれらと似たような行がある はずなので、これもヒントにして適切な入れ場所を考えてほしい。 もしかしたら既に気付いている人もいるだろうが、各行の先頭についている「スペース」の数もヒントとして残しているつもりである。

x の平均と分散を求める

必要な手順 (の要点だけ) を以下に示す。

• 「西暦の分散」を計算するのに必要な器を用意する。 例えば、「総和」用に seireki_souwa、「2乗の総和」用に seireki2_souwa、さらに「平均」、「2乗の平均」「分散」用に seireki_heikin、seireki2_heikin、seireki_bunsan をそれぞれ宣言する。 ただし、少なくとも 「平均」、「2乗の平均」「分散」を入れる変数は実数型変数として宣言 しておくこと。 その理由は、もとになる「西暦」は整数であるとはいえ、割り算したり平方根をとったりしたら整数でなくなる可能性がある からである。

  i

  n

  t

  e

  g

  e

  r

  :

  :

  s

  e

  i

  r

  e

  k

  i

  _

  s

  o

  u

  w

  a

  ,

  s

  e

  i

  r

  e

  k

  i

  2

  _

  s

  o

  u

  w

  a

  r

  e

  a

  l

  :

  :

  s

  e

  i

  r

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

  ,

  s

  e

  i

  r

  e

  k

  i

  2

  _

  h

  e

  i

  k

  i

  n

  ,

  s

  e

  i

  r

  e

  k

  i

  _

  b

  u

  n

  s

  a

  n

• 西暦の「総和」と「2乗の総和」を入れる「器」の中身をあらかじめ「空っぽ」にしておく。

  s

  e

  i

  r

  e

  k

  i

  _

  s

  o

  u

  w

  a

  =

  0

  s

  e

  i

  r

  e

  k

  i

  2

  _

  s

  o

  u

  w

  a

  =

  0

• 西暦の「総和」と「2乗の総和」を入れる「器」のに値をどんどん足していく。

  d

  o

  r

  e

  a

  d

  (

  *

  ,

  '

  (

  i

  4

  ,

  1

  2

  f

  5

  .

  1

  )

  '

  ,

  i

  o

  s

  t

  a

  t

  =

  j

  o

  u

  t

  a

  i

  )

  s

  e

  i

  r

  e

  k

  i

  ,

  k

  i

  o

  n

  i

  f

  (

  j

  o

  u

  t

  a

  i

  /

  =

  0

  )

  e

  x

  i

  t

  s

  e

  i

  r

  e

  k

  i

  _

  s

  o

  u

  w

  a

  =

  s

  e

  i

  r

  e

  k

  i

  _

  s

  o

  u

  w

  a

  +

  s

  e

  i

  r

  e

  k

  i

  s

  e

  i

  r

  e

  k

  i

  2

  _

  s

  o

  u

  w

  a

  =

  s

  e

  i

  r

  e

  k

  i

  2

  _

  s

  o

  u

  w

  a

  +

  s

  e

  i

  r

  e

  k

  i

  *

  *

  2

  e

  n

  d

  d

  o

• 西暦の「平均」と「分散」を計算する。 ただし 前回 にも示したように、seireki_souwa や nensuu のような整数型変数の値を実数だと思って計算してほしいところでは real(実数でない変数) とする こと。

  s

  e

  i

  r

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

  =

  r

  e

  a

  l

  (

  s

  e

  i

  r

  e

  k

  i

  _

  s

  o

  u

  w

  a

  )

  /

  r

  e

  a

  l

  (

  n

  e

  n

  s

  u

  u

  )

  s

  e

  i

  r

  e

  k

  i

  2

  _

  h

  e

  i

  k

  i

  n

  =

  r

  e

  a

  l

  (

  s

  e

  i

  r

  e

  k

  i

  2

  _

  s

  o

  u

  w

  a

  )

  /

  r

  e

  a

  l

  (

  n

  e

  n

  s

  u

  u

  )

  s

  e

  i

  r

  e

  k

  i

  _

  b

  u

  n

  s

  a

  n

  =

  s

  e

  i

  r

  e

  k

  i

  2

  _

  h

  e

  i

  k

  i

  n

  -

  s

  e

  i

  r

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

  *

  *

  2

まずここまで作ったところで再び、gfortran にかけ、エラーが出ないことを確かめよ。

また念のため、プログラムの最後のほうで

w

r

i

t

e

(

*

,

*

)

s

e

i

r

e

k

i

_

h

e

i

k

i

n

,

s

e

i

r

e

k

i

2

_

h

e

i

k

i

n

,

s

e

i

r

e

k

i

_

b

u

n

s

a

n

のような記述を追加した上で gfortran にかけ、プログラムを実行させてみることにより、西暦の「平均」や「分散」が正しく求まっているかを確認してみるのもよい「お作法」である。 ちなみに、西暦の「平均」は \((1890+2023)/2=1956.500000\)、「分散」は \((134^2-1)/12=1496.250000\) みたいになるのが正解 (手計算でやってみるのも吉)。 なお、西暦の「平均」や「分散」が正しく求まっていることが確認できたら、この write 文は (もう使わないから) 消してしまっても構わない。

x と y の積の平均を求める

ここでも必要な手順 (の要点だけ) を以下に示す。

• 「気温と西暦の積の平均」を求めるのに必要な「器」を用意する。 例えば「総和」用に seki_souwa、「平均」用に seki_heikin という名前の 実数 を宣言する。

  r

  e

  a

  l

  :

  :

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  ,

  s

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

• 「気温と西暦の積の総和」を入れる「器」の中身をあらかじめ「空っぽ」にしておく。

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  =

  0

  .

  0

• 「気温と西暦の積の総和」を入れる「器」にどんどん値を足していく。 各年の各月の気温を読むたびに、seki_souwa という器に足し込んでいる。 ただし、平均気温 (実数) と西暦 (整数) のかけ算をする際には安全のため、seireki* ではなく real(seireki)* とする こと。

  d

  o

  r

  e

  a

  d

  (

  *

  ,

  '

  (

  i

  4

  ,

  1

  2

  f

  5

  .

  1

  )

  '

  ,

  i

  o

  s

  t

  a

  t

  =

  j

  o

  u

  t

  a

  i

  )

  s

  e

  i

  r

  e

  k

  i

  ,

  k

  i

  o

  n

  i

  f

  (

  j

  o

  u

  t

  a

  i

  /

  =

  0

  )

  e

  x

  i

  t

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  =

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  +

  r

  e

  a

  l

  (

  s

  e

  i

  r

  e

  k

  i

  )

  *

  k

  i

  o

  n

  (

  5

  )

  e

  n

  d

  d

  o

• 「気温と西暦の積の平均」を計算する。 ここでも当然、/nensuu ではなく /real(nensuu) とする 必要がある。

  s

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

  =

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  /

  r

  e

  a

  l

  (

  n

  e

  n

  s

  u

  u

  )

次にここまで作ったところで、gfortran にかけ、エラーが出ないことを確かめよ。

「あてはめ」られた直線の方程式を求める

ここまでくれば後は簡単、のはず。 必要な手順 (の要点だけ) を示す。

• 必要な「器」を用意する。 例えば 直線の傾きを kaiki_a、y切片を kaiki_b に入れることにする。

  r

  e

  a

  l

  :

  :

  k

  a

  i

  k

  i

  _

  a

  ,

  k

  a

  i

  k

  i

  _

  b

• 直線の傾き \(a\) と y切片 \(b\) の値を計算し、表示させる。 講義のノートを参照して、適切な計算式を書き入れよ。 その際、何が x (横軸) で、何が y (縦軸) か に十分注意すること。 (わざわざ文字にして書いてみると当たり前すぎるけれど) 西暦の平均 seireki_heikin と5月の気温の平均 kion_tsuki_heikin(5) は違うものだし、同様に西暦の分散 seireki_bunsan と5月の気温の分散 kion_tsuki_bunsan(5) も違うものなのだから。 また \(a\) と \(b\) の値を出力するにあたり、傾き \(a\) は (値が小さいから) 小数点以下を多め (この例では小数点以下4桁) に表示するようにしてある。 詳しい説明は ここ を参照のこと。

  k

  a

  i

  k

  i

  _

  a

  =

  k

  a

  i

  k

  i

  _

  b

  =

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  f

  8

  .

  4

  ,

  f

  6

  .

  1

  )

  '

  )

  k

  a

  i

  k

  i

  _

  a

  ,

  k

  a

  i

  k

  i

  _

  b

ここまでを参考にして、5月の平均気温の過去134年間の変化を「あてはめ」た直線の方程式を求める手順がどう書けるかを考えよ。 そして前回に作成したプログラムにその機能を追加したプログラムを作成せよ [みほん]。

gfortran に怒られないプログラムができたら、例の如く

$ cat matsuyama-kion.txt | ./a.out

と実行して計算結果を確認せよ。 正しくプログラミングできていたら、プログラムの出力 (の最後の1行) は以下のようになっているはずである。

  0.0246 -30.1

最小二乗法による直線への「あてはめ」 (その2)

回帰曲線の方程式が求まったら、それを実際にグラフに表示して確かめてみよう。 端末から Gnuplot を起動し、以下のように入力してみよう (gnuplot> は Gnuplot のプロンプト)。

gnuplot> plot "matsuyama-kion.txt" using 1:6,0.0246*x-30.1

とすると、「5月の平均気温」の年変化を表わすグラフとともに、いま求めた回帰直線が緑で示されているはずだ。 直線への「回帰」がどの程度うまくいっているかを確認してみよ。 さらに、

gnuplot> plot "matsuyama-kion.txt" using 1:6,0.0246*x-30.1,17.97

として、前回までに求めた5月の平均気温も合わせて表示してみよ。 その際、2つの直線の位置関係はどうなっているか。

決定係数を求める

最小二乗法による「あてはめ」の程度を調べる最初の例として、5月の平均気温の過去134年間の変化の「あてはめ」の決定係数 \(r^2\) を求めてみよう。 必要な手順 (の要点だけ) を以下に示す。 ただし、すでに「あてはめ」られた直線の方程式が求まっている ので、新たに計算しなければいけないものは何もない。

• 「決定係数」を求めるのに必要な「器」を用意する。 例えば r2 という名前の実数を宣言する。

  r

  e

  a

  l

  :

  :

  r

  2

• 決定係数 \(r^2\) を計算し、表示させる。 講義のノートを参照して、適切な計算式を書き入れよ。 ただし5月の平均気温の回帰直線の傾き \(a\) (この例では kaiki_a) が既に求まっているので、kaiki_a を使った式で r2 を計算するのが簡単 でよい。 また r2 の値を出力するにあたり、(傾き \(a\) と同様に値が小さいから) 小数点以下を多め (この例では小数点以下4桁) に表示するようにしてある。 詳しい説明は ここ を参照のこと。

  		r	2		=
  		w	r	i	t	e	(	*	,	'	(	f	8	.	4	,	f	6	.	1	,	f	7	.	3	)	'	)	k	a	i	k	i	_	a	,	k	a	i	k	i	_	b	,	r	2

これらを参考にして、5月の平均気温の過去134年間の変化を直線で「あてはめ」た際の決定係数を求める手順がどう書けるかを考えよ。 そして上で作成したプログラムをもとに、決定係数を計算する機能を追加したプログラムを作成せよ [みほん]。

gfortran に怒られないプログラムができたら、例の如く

$ cat matsuyama-kion.txt | ./a.out

と実行して計算結果を確認せよ。 正しくプログラミングできていたら、プログラムの出力 (の最後の4行) は以下のようになっているはずである。

 134  5.26  5.59  8.64 13.72 17.97 21.91 26.24 27.05 23.48 17.69 12.51  7.68
 134  1.51  2.10  2.00  1.52  1.33  0.98  1.37  1.09  1.74  1.53  1.60  1.53
 134  1.23  1.45  1.42  1.23  1.15  0.99  1.17  1.04  1.32  1.24  1.27  1.24
  0.0246 -30.1  0.679

最小二乗法による直線への「あてはめ」 (その3)

回帰曲線の方程式とその決定係数が求まったら、「あてはめ」がどの程度うまくいっているのかを実際にグラフに表示して確かめてみよう。 端末から Gnuplot を起動し、以下のように入力してみよう (gnuplot> は Gnuplot のプロンプト)。

gnuplot> plot "matsuyama-kion.txt" using 1:6,0.0246*x-30.1,0.0246*x-30.1+sqrt(1-0.679)*1.15,0.0246*x-30.1-sqrt(1-0.679)*1.15

これによって得られたグラフと、前回 で

gnuplot> plot "matsuyama-kion.txt" using 1:6,17.97,17.97+1.15,17.97-1.15

とした際に得られたグラフとを比較して、最小二乗法の決定係数の意味を確認せよ。 例えば (非常に長いので、「コピペ」するのがよい)

gnuplot> plot "matsuyama-kion.txt" using 1:6,0.0246*x-30.1,0.0246*x-30.1+sqrt(1-0.679)*1.15,0.0246*x-30.1-sqrt(1-0.679)*1.15,17.97,17.97+1.15,17.97-1.15

として表示される「ななめの『帯』」と「水平の『帯』」の2つの帯の幅の違いは何を意味しているか。

最小二乗法による直線への「あてはめ」 (その4)

今度は、全ての月の平均気温の年変化をあてはめた直線の方程式と、あてはめの決定係数をいっぺんに求めるプログラム enshu04c.f90 を作ってみよう [みほん]。 なお「完成までの途中経過」にあたるものとして、(12ヶ月分を求める準備はしてあるのに) 5月の気温だけを最小二乗法するプログラム例を ここ に載せておく。

必要な手順 (の要点の一部) を以下に示す。

• 必要な「器」を定義する。 「気温と西暦の積の平均」を求めるのに必要な変数、および直線の傾きと y 切片や決定係数にあたるものを、12個の実数を入れる配列 として宣言し直す。

  r

  e

  a

  l

  ,

  d

  i

  m

  e

  n

  s

  i

  o

  n

  (

  1

  2

  )

  :

  :

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  ,

  s

  e

  k

  i

  _

  h

  e

  i

  k

  i

  n

  r

  e

  a

  l

  ,

  d

  i

  m

  e

  n

  s

  i

  o

  n

  (

  1

  2

  )

  :

  :

  k

  a

  i

  k

  i

  _

  a

  ,

  k

  a

  i

  k

  i

  _

  b

  ,

  r

  2

• 各月の「平均気温と西暦の積の平均」を計算し、さらに「あてはめ」た直線の傾きと y 切片や決定係数を求める。 しつこいようだが、まさかこの期に及んで

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  )

  =

  0

  .

  0

  から

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  2

  )

  =

  0

  .

  0

  まで、あるいは

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  )

  =

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  )

  +

  r

  e

  a

  l

  (

  s

  e

  i

  r

  e

  k

  i

  )

  *

  k

  i

  o

  n

  (

  1

  )

  から

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  2

  )

  =

  s

  e

  k

  i

  _

  s

  o

  u

  w

  a

  (

  1

  2

  )

  +

  r

  e

  a

  l

  (

  s

  e

  i

  r

  e

  k

  i

  )

  *

  k

  i

  o

  n

  (

  1

  2

  )

  まで、さらには

  k

  a

  i

  k

  i

  _

  a

  (

  1

  )

  =

  k

  a

  i

  k

  i

  _

  b

  (

  1

  )

  =

  r

  2

  (

  1

  )

  =

  から

  		k	a	i	k	i	_	a	(	1	2	)		=
  		k	a	i	k	i	_	b	(	1	2	)		=
  		r	2	(						1	2	)		=

  の如く、12個の式を長々と書き並べるなんてことのないように。

• 最後に結果を表示する。 この部分は以下のようにするのがよい。 これは \(a\) の値 (正確にはその絶対値) が小さいので、10倍した値を小数点以下3桁まで表示する のがちょうどいい感じになるからである。 また 4x とは「空白を4つ入れてください」という意味。 詳しい説明は ここ を参照のこと。

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  i

  4

  ,

  1

  2

  f

  6

  .

  2

  )

  '

  )

  n

  e

  n

  s

  u

  u

  ,

  k

  i

  o

  n

  _

  t

  s

  u

  k

  i

  _

  h

  e

  i

  k

  i

  n

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  i

  4

  ,

  1

  2

  f

  6

  .

  2

  )

  '

  )

  n

  e

  n

  s

  u

  u

  ,

  k

  i

  o

  n

  _

  t

  s

  u

  k

  i

  _

  b

  u

  n

  s

  a

  n

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  i

  4

  ,

  1

  2

  f

  6

  .

  2

  )

  '

  )

  n

  e

  n

  s

  u

  u

  ,

  k

  i

  o

  n

  _

  t

  s

  u

  k

  i

  _

  h

  e

  n

  s

  a

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  4

  x

  ,

  1

  2

  f

  6

  .

  3

  )

  '

  )

  k

  a

  i

  k

  i

  _

  a

  *

  1

  0

  .

  0

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  4

  x

  ,

  1

  2

  f

  6

  .

  1

  )

  '

  )

  k

  a

  i

  k

  i

  _

  b

  w

  r

  i

  t

  e

  (

  *

  ,

  '

  (

  4

  x

  ,

  1

  2

  f

  6

  .

  3

  )

  '

  )

  r

  2

正しくプログラミングできていれば、プログラムの出力 (の最後の6行) は以下のようになるはずである。 先と同様のグラフを11月以外の月 (例えば1月) の平均気温についても描いてみることにより、決定係数の大きい場合と小さい場合とで「あてはめ」の程度がどのように異なっているか確認せよ。

 134  5.26  5.59  8.64 13.72 17.97 21.91 26.24 27.05 23.48 17.69 12.51  7.68 <- 平均
 134  1.51  2.10  2.00  1.52  1.33  0.98  1.37  1.09  1.74  1.53  1.60  1.53 <- 分散
 134  1.23  1.45  1.42  1.23  1.15  0.99  1.17  1.04  1.32  1.24  1.27  1.24 <- 標準偏差
     0.146 0.200 0.230 0.198 0.246 0.156 0.154 0.181 0.195 0.240 0.202 0.141 <- 傾き a の10倍
     -23.3 -33.5 -36.4 -25.1 -30.1  -8.6  -4.0  -8.3 -14.7 -29.2 -27.1 -19.9 <- y 切片 b
     0.211 0.285 0.395 0.387 0.679 0.373 0.260 0.449 0.328 0.562 0.383 0.196 <- 決定係数 r2

正しく動いていることが確認できたら、本日の演習の出席確認のため、このプログラムを 末尾のフォーム 経由で提出せよ。

出席確認用プログラム提出フォーム

Moodle のコースページ 以下の 課題ファイル提出 により、各自で作成した (enshu04c.f90 に相当する) 過去134年間の各月の気温の変化を最小二乗法で直線にあてはめる Fortran 90 プログラムの最終版を提出せよ。