02. グラフを描こう

Photo by Chris Liverani on Unsplash

計算結果は、数字を眺めるだけでは特徴を見落としやすい。今日は Julia の Plots package を使い、何を横軸・縦軸に置いたか説明できる図を作ろう。コマンドの暗記よりも、「予想する → 描く → 読み取る」を繰り返すことを重視する。

本時の道筋（90分）

この授業で試す作図方法の全体像を先に見ておこう。関数・ベクトル・行列のどれを渡すかを意識しながら、各節の例に進む。

授業前の準備。 Julia/Jupyter を起動し、まだ入れていなければ次を一度だけ実行しておこう。Plots は大きな package なので、初回のインストールや最初の描画には時間がかかることがある。

1using Pkg
2Pkg.add("Plots")

インストール済みか確認したければ using Pkg; Pkg.status() を実行する。同じ package を再度指定してもよいが、授業時間を使わないよう、前もって確認するのがよい。利用するセッションでは、最初に次を実行する。

1using Plots

1. 関数を描く：まず形を予想する

plot は関数を渡すと、その値をいくつかの点で評価して曲線を描く。まず sin の山と谷がどう現れるか予想してから実行しよう。

1plot(sin)

描画範囲も自分で決められる。次の -π から π には、零点はいくつあるだろうか。

1plot(sin, -π, π; xlabel="x", ylabel="sin(x)", label="sin(x)")

名前のない関数も、x -> ... と書けばその場で渡せる。元の例の cos(x) * x^2 は、左右対称だろうか。

1plot(x -> cos(x) * x^2)

試す（2分） x^2 を x に替えて再描画する。対称性はどう変わったか、式と図を結びつけて一文で説明しよう。

2. データを描く：横軸は何か

次は関数ではなく、計算で得た数の並びを描こう。元の例と同じく、$-2$ から $2$ まで 0.1 刻みで二乗した値を作る。

1v = [n^2 for n in -2:0.1:2]
2length(v), eltype(v)

結果は「41 個の Float64」のはずだ。表示される各数値を全部読む必要はない。要素数・型・最小値・最大値に注目しよう。

1minimum(v), maximum(v)
2plot(v)

ここでの横軸は -2:0.1:2 ではなく、要素番号 1, 2, …, 41 である。元の $n$ を横軸にするなら両方を渡す。

1x = collect(-2:0.1:2)
2plot(x, v; xlabel="n", ylabel="n²", label="n²")

試す（4分） plot(v) と plot(x, v) で、谷の位置を横軸の値として読み取ろう。同じ縦軸データでも、横軸の意味を省くと解釈が変わる。

測定点を見せたい場合はマーカーを付ける。marker = :auto と marker = :hexagon を比べよう。

1plot(v, marker=:auto)

1plot(v, marker=:hexagon)

Julia でベクトルは 1 次元配列、行列は 2 次元配列である。plot に行列だけを渡すと、列ごとに別の系列として描かれる。この性質は次の節で確かめる。

3. 行列を描く：何を色や高さで表すか

同じデータでも、等高線、色、立体表示では読み取りやすい性質が異なる。ここでは $z=x^2+y^2$ の値を格子点に並べる。

1A = [x^2 + y^2 for x in -1:0.1:1, y in -1:0.2:1]
2size(A), minimum(A), maximum(A)

size(A) は (21, 11)。長い行列出力を全部読む代わりに、中央付近が小さく、端に向かって大きくなることを予想してみよう。以下の図は**すべて同じ A**から作ったものだ。冒頭の全体図では「行列」からどの表示方法につながっていたか、思い出そう。

高さを線・色で見る

contour(A) は同じ高さの場所を線で結ぶ。どこが最小で、中心から離れるとどう変わるかを図から読もう。

1contour(A)

一方、plot(A) は行列の各列を別の系列として扱う。この場合は 11 本の曲線になる。

1plot(A)

heatmap(A) では高さを色に置き換える。contour との違いは何だろうか。

1heatmap(A)

等高線の間を塗るなら fill = true とする。

1contour(A, fill=true)

試す（5分） contour(A) と heatmap(A) を実行し、最小値の場所を両方で指し示そう。どちらが「同じ高さの領域」を読みやすいだろうか。

表示を重ねる場合、! が付くコマンドは今の図に要素を足す。色と線を重ねる次の例も試せる。

1heatmap(A, opacity=0.85)
2contour!(A)

高さを立体で見る

wireframe は格子の形を、surface は面の形を見せる。Web ブラウザや PC によって表示に時間がかかるため、授業中はまず図を見比べ、余裕があれば自分でも実行しよう。

1wireframe(A)

1surface(A)

3 次元表示の見る方向は camera=(平面角, 高さ角) で変えられる。元の (22, 58) の例では、山や谷の見え方がどう変わるだろうか。

1wireframe(A, camera=(22, 58))

4. 人に見せる図にして保存する

軸名・凡例・図の大きさを与えてから保存しよう。ここでは 4 次関数を使う。

1p = plot(x -> x^4 - 20x^2 + 1;
2         xlabel="x", ylabel="f(x)", label="f(x)",
3         size=(1200, 800))
4savefig(p, "sample.png")

もとのように、描いた後から plot!(size=(1200,800)) で大きさを変え、savefig("sample.png") として現在の図を保存することもできる。sample.png は Julia の現在の作業フォルダーにできる。Jupyter のファイラーでも確認しよう。

試す（3分） sample.png を開き、軸名・凡例・曲線が読めるか確認する。「画面に出た」だけで終わらせず、保存したファイルまで確かめよう。

保存形式は backend に依存する。標準設定で用いられる GR では、PNG のほか PDF、PS、SVG などが使える。画面とファイルの関係を次の図でも確認しておこう。

5. 短い GIF を試す

@gif はループ内で作った図をフレームとして集め、GIF にする。まずは授業中に試しやすいよう、元の 91 フレームの例を 10 フレーム に縮めてみよう。初回は動画用の外部ツールの準備に時間がかかる場合がある。

1@gif for t in 0:10:90
2    surface(A, camera=(t, 30))
3end

実行後に表示される GIF または保存先のメッセージを確認しよう。@gif の保存先は実行環境によって異なる。以下は 91 フレーム版で、カメラ角を少しずつ変えている。

1@gif for t in 0:1:90
2   surface(A, camera=(t,30))
3end

試す（3分） 10 フレーム版ができたら、0:10:90 を 0:15:90 に変える。フレーム数と動きの滑らかさ・作成時間の関係を観察しよう。

6. 今日の小課題（残り8分）

次のデータは $e^{-x^2}$ を格子点で計算したものである。上の例を参考に、横軸を x、縦軸を計算値として図を作り、軸名と凡例を付けて PNG に保存しよう。

1x = range(-2, 2; length=81)
2y = exp.(-x.^2)

できたら、(1) 最大値はどこか、(2) 左右対称か、(3) なぜ plot(y) だけでは横軸の説明として不十分か、を一文ずつ書こう。授業終了時には、図を見て「どの変数をどの軸に置いたか」説明できればよい。

レポート No.02（授業後に取り組む）

課題 1・2 では授業中の小課題と異なるデータを自分で選び、図が何を示すか説明すること。

1. 授業用 web に出てくる例と異なる 2 次元グラフ、3 次元グラフを作り、データの解説とグラフを提出する。
2. 授業用 web に出てくる例と異なる GIF グラフを作り、レポートに添付して提出する。
3. 標本点の間を見落とさない。 $f(x)=\sin(5x)$ を $0\leq x\leq 2\pi$ で描く。まず 11 点だけを使った図を作り、次に 201 点で描いた曲線と 11 点を同じ座標軸に重ねる。2 枚の図を示し、最初の図からどのような誤解が生じるか、標本点の数を増やすと何が分かるかを説明する。凡例で 2 種類の点を区別する。
4. 刻み幅を小さくすれば常によいか。 $f(x)=\sin x$ の $x=1$ における微分を前進差分 $D(h)=[f(1+h)-f(1)]/h$ で近似する。$h=10^{-1},10^{-2},\ldots,10^{-12}$ に対する誤差 $|D(h)-\cos 1|$ を両対数グラフにし、誤差が最小となる刻み幅を読み取る。小さい $h$ で誤差が再び大きくなる理由を考察する。横軸を刻み幅 $h$、縦軸を誤差として軸名を付ける。

課題 3 の標本点は x11 = range(0, 2π; length=11)、x201 = range(0, 2π; length=201) で作れる。関数を配列に適用するときは sin.(5 .* x11) のように書く。課題 4 では h = 10.0 .^ (-1:-1:-12) として刻み幅を作り、前進差分を計算してから abs. で誤差を求める。図だけでなく、読み取ったことを自分の言葉で添えること。

レポートは 学籍番号-氏名-02.pdf とし、web フォームから提出する。実行形式のプログラムをメールに添付すると受信側に拒否される場合があるため、PDF を使う。紙媒体の手渡しも可能だが、破損・紛失の可能性がある。

もっと調べたいとき

授業ではコマンドの全オプションを覚える必要はない。困ったときに調べる入口を覚えよう。

1?plot
2plotattr(:Series)
3plotattr("markershape")
4plotattr("camera")

?plot は plot と plot! の違い、plotattr(:Series) は系列の属性、plotattr("markershape") はマーカーの選択肢を教えてくれる。詳細は Plots のチュートリアル と Animations を参照しよう。