表計算ソフトの基礎 31 5 最小二 ( 自 ) 乗法 (Least squares method) 理工系の実験などでは, たいてい測定が行われる. ただし, 測定コストを減らすため, 粗い測定を行い, 得られた ( 誤差を含む ) 測定点群を説明することができる, 真の特性として何らかの関数 (

Transcription

1 表計算ソフトの基礎 31 5 最小二 ( 自 ) 乗法 (Least squares method) 理工系の実験などでは, たいてい測定が行われる. ただし, 測定コストを減らすため, 粗い測定を行い, 得られた ( 誤差を含む ) 測定点群を説明することができる, 真の特性として何らかの関数 ( モデル ) を想定し, 測定点群に近似曲線を当てはめる ( フィッティング ) という方法がよく用いられる. ここで, 実験により,N 個の測定点群 D = {d 1,,d N } が得られたとする. 具体的には, 例えば, 二次元のデータ d k = (x k, y k ) が N=8 個得られて, 次のような値となっているとする. (x 1, y 1 )=(0, 40), (x 2, y 2 )=(5, 25), (x 3, y 3 )=(6, 70), (x 4, y 4 )=(11, 55), (x 5, y 5 )=(12, 100), (x 6, y 6 )=(17, 85), (x 7, y 7 )=(18, 130), (x 8, y 8 )=(23, 115) これらのデータ D を, 表計算ソフトを用いてグラフ化すると, 例えば, 図 61 のようになる. これらの測定点群は, 誤差を含んでいるとすれば, 当てはめる曲線が必ずしも測定点を通る必要はない. このような状況において, 表計算ソフトには, 何らかのデータ群に対して近似曲線を自動的に当てはめる機能が用意されている場合がある.Calc 図 61 データ D をプロットした例や Excel にもその機能がある. ここでは, その機能を試してみることにしよう. 5.1 表計算ソフトを用いて近似曲線を当てはめる方法 データのプロットまず, データ D を, 図 62 左側のように入力してみよう. ただし,A 列に記載した (x k, y k ) の表記は記載しなくても良い. そして,B2:C10 をアクティブ化して, 右のようなアイコンをクリック, もしくはメニューの 挿入 (I) グラフ(A) を選択すると, 図 62 右のような グラフウィザード が現れる. ここで, 散布図 の 点のみ を選び, 完了 ボタンを押した結果が図 61 である. 図 62 データの入力例 ( 左 ) とグラフウィザードウィンドウ ( 右 )

2 表計算ソフトの基礎 近似曲線の挿入 では, 図 61 のグラフに近似曲線を当てはめる機能を試してみることにする. まず 1. グラフをダブルクリックして編集可能状態にする. 2. 任意のデータ点上にマウスのカーソルを移動させ, 左クリックして選択状態にする. 3. データ点が選択状態にあるとき, マウスを右クリックすると図 63 左のようになる. 4. サブメニューの中から図 63 右のように 近似曲線を挿入 を選択する. 図 63 データ点群に近似曲線を当てはめるサブウィンドウの操作 すると, 図 64 のような近似曲線を当てはめるための専用ウィンドウが現れる. 選択できる近似曲線は, 以下の通りである. 線形 : 直線を当てはめる. 対数 : 対数関数を当てはめる. 指数 : 指数関数を当てはめる. 累乗 : 累乗関数を当てはめる. 多項式 : 多項式を当てはめる. 移動平均 : 移動平均を当てはめる. 近似曲線の中に, 直線も含まれていることがわかる. データ D に, 近似曲線を当てはめた例を図 65 に示す. 図 65 には, 図 64 近似曲線当てはめ専用ウィンドウ (a) 対数関数 (b) 累乗関数 (c) 指数関数 (d) 多項式を当てはめた結果が示されている. ただし, 多項式については, 次数 という設定項目があるが, 図 65(d) では, 次数に 3 を設定した結果を示している. どの曲線 ( モデル ) が尤もらしいかについては, 曲線を当てはめる者の主観が入り得るため, データを扱う者のセンスが問われる. 図 65 近似曲線を当てはめた例

3 5.1.3 挿入した近似曲線の削除挿入した近似曲線を削除したときは, 以下のようにすればよい. 1. 挿入された近似曲線上をマウスで左クリック 図 66 のように, 近似曲線が選択状態になる. 2. 近似曲線が選択された状態でマウスを右クリック 図 66 のように, サブメニューが出る 近似曲線を削除(N) をクリック 表計算ソフトの基礎 33 図 66 削除する近似曲線の選択とサブメニューの選択 練習 に示したいずれかの近似曲線を挿入し, 削除してみなさい. 5.2 直線フィッティング 5.1 で示したように, データに近似曲線を当てはめる機能が簡単に使えることがわかる. 一般的に, このようなデータに対する曲線フィッティングにおいて, 最も基本的な方法と言えるのが最小二乗法 (Least squares method) であり, 最小自乗法とも呼ばれる. ここでは, 最小二乗法の中でも, 最も単純な, 直線フィッティングを通じて最小二乗法を概観する 表計算ソフトを用いて最も簡単に直線を当てはめる方法表計算ソフトを用いて測定点群に直線を当てはめる最も簡単な方法は, で示した近似曲線の挿入における, 図 64 のウィンドウにおいて, 線形 を選ぶことである. すると, 図 67 のような結果が得られる. この線形近似直線はどのような関数であるかを知りたい場合, その関数 f(x) を表示させる機能も存在する. 関数 f(x) を表示させるには, 挿入された線形近似直線上にマウスポインタを図 67 直線を当てはめた例

4 表計算ソフトの基礎 34 図 68 直線を当てはめた例 置き, 次のようにすればよい. 1. 挿入された線形近似直線上をマウスで左クリック 図 68 左のように, 線形近似直線が選択状態になる. 2. 線形近似直線が選択された状態でマウスを右クリック 図 68 左のように, サブメニューが出る 近似曲線の方程式を挿入(I) をクリック すると, 図 68 右の赤枠のように, f(x) = x という関数表記が現れる. これは, 傾き (slope) a : 切片 (intercept) b : とする直線の方程式 y = f (x) = a x + b を表している. 実は, データ D からこの線形近似直線を求めるため, 傾き a と切片 b を算出する方法と して用いられているのが最小二乗法である 最小二乗法で傾きと切片を算出する表計算ソフトの関数 1 Calc( や Excel) には, 近似直線 ( 回帰直線 ) の傾きと切片をそれぞれ別々に計算する関数がある. a = SLOPE(Data_Y, Data_X) b = INTERCEPT(Data_Y, Data_X) 例えば, 図 62 のシートでは,N=8 個のデータ d k = (x k, y k ) の成分 x k が B3:B10,y k が C3:C10 に格納されていた. このとき, 例えば空いているセル E3 に =SLOPE(C3:C10,B3:B10) と入力すれば, そのセル値が と表示される. また, 切片も, 例えばセル E4 =INTERCEPT(Data_Y, Data_X) と入力すれば, そのセル値が と表示される. その様子を示したものが図 69 である の傾き 切片の値に一致している. 図 69 傾きと切片の算出結果

5 表計算ソフトの基礎 最小二乗法で傾きと切片を算出する表計算ソフトの関数 2 <LINEST 関数と INDEX 関数を使う方法 > Calc( や Excel) には, 直線の傾きおよび切片を計算できる関数として LINEST がある. 計算結果 = LINEST(Range1 (, Range2) (,Option1) (,Option2)) Range1 ( 必須 ) : 従属関数 (y) の行列 Range2 ( 任意 ) : 独立関数 (x) の行列 Option1 ( 任意 ) : y 切片を 0 に指定するか否かの指定 ( 省略時は True) 1 (True) y 切片を計算する 0 (False) y 切片を 0 とする Option2 ( 任意 ) : 回帰直線の補正項を返すか否かの指定 ( 省略時は True) 1 (True) 回帰直線の補正項は返す 0 (False) 回帰直線の補正項は返さない Range1 と Range2 のセル数が一致していない場合はエラーとなる 例えば,5.2.2と同じような結果を得るには, 例えばセル E6 に =LINEST(C3:C10, B3:B10) とすれば,Option1 としては y 切片を計算で求め, 回帰直線の補正項は返さないので, が得られる. これは,5.2.1の傾きの値と一致する. では, 切片はどうなっているだろうか. 実は,LINEST が返す計算結果は, 単一の数値ではなく, 構造を持っている. =LINEST(C3:C10, B3:B10) のようにセルへ式を入力している場合は, } のように, 傾きの値と切片の値が,2 行 1 列で行列表現されているのであるが, セルの表示では, 先頭 の要素 (1 行 1 列の要素 ) のみが表示されるようになっている. このような行列表現が返されるとき, 各要素の添え字 (Index) を指定して取り出す関数が用意されている. 要素の値 = INDEX( 参照 (, 行 ) (, 列 ) (, 範囲 )) 参照 ( 必須 ) : 参照の対象 ( 行列表現など ) や範囲名 行 ( 任意 ) : 何行目かを指定. 指定がない (0) の場合はすべての行が返る 列 ( 任意 ) : 何列目かを指定. 指定がない (0) の場合はすべての列が返る 範囲 ( 任意 ) : 多重範囲を参照する場合の部分範囲の Index を表す 今回の事例では, 列が一つであるため, 列の指定は必要でなく, 範囲も無視すれば, 例えばセル E7 に =INDEX(LINEST(C3:C10, B3:B10),2) と 2 行目を指定することで切片の値を取り出すことができる. ちなみに, =INDEX(LINEST(C3:C10, B3:B10),1) は, =LINEST(C3:C10, B3:B10) と入力することと同じことになる. 練習 に示した線形近似直線を実際に表示させてみよう. また,5.2.2 の SLOPE および INTERCEPT 関数, また,5.2.3 に示した LINEST および INDEX 関数を用いて線形近似直線の傾きと切片を求めてみよう.

6 表計算ソフトの基礎 36 <LINEST 関数を用い, 傾きと切片を同時に得る方法 > 実は,LINEST が返してくる行列表現に対し,INDEX を使わず, 傾き 切片を同時に取り出す方法がある. 適当な連続した二つのセル, 例えば E9:F9 をアクティブ化し, =LINEST(C3:C10, B3:B10) を数式入力枠に入力する. そして, Ctrl キー + Shift キー + Enter キーを押せばよい. すると, 図 70 のようになる. 図 70 の E9 には傾き,F9 には切片の値が同時に得られていることがわかる. このとき, E9:F9 は, 連続した領域として特殊な状態にある. 図 70 では F9 が選択されているが, その中身は図 70 の数式入力枠に見えている {=LINEST(C3:C10,B3:B10)} のように 式が括弧 { } で囲まれていることがわかるだろうか. これは配列数式と呼ばれるものである.E9 を選択状態にしても, 数式図 70 Ctrl + Shift + Enter で式を入力した結果入力枠には同じ式が現れる. つまり, 配列数式は E9:F9 の連続したセルに対応づけられている. ここでは詳しく述べないが, 表計算ソフト上で, 行列演算を行うときに良く使われる手法である 傾きおよび切片から線を描く 5.2.1で示したように, データから直線をグラフに描画する機能が表計算ソフトに用意されていることは示したが, 線形近似直線の傾きおよび切片が得られたとき, 線形近似直線を元のデータと一緒のグラフに描くにはどうすればよいのだろうか. それは, y = a x + b の右辺において,a,b が定まっているので, 線を引こうとする x の範囲において,2 点以上, 各 x に対応する y を表計算上で算出し, 元データと一緒にグラフ化したのち,(x, y ) の系列について, 線形近似直線を挿入すればよい.y の系列は, もともと直線上の点系列であるから, 線形近似直線は, データ点 (x, y ) を通る直線となる. 例えば, 図 70 のセル B3 を x とする y を求めるため,D3 に, = * B と入力する. ただし, この入力法では, いちいち a,b の値を入力する必要がある. しかし, 図 70 の E3 に a,e4 に b の値があるため, これを常に使うならば, = $E$3 * B3 + $E$4 とすればオートフィルを考慮した式となる. 直線を引くには2 点必要であるため, この D3 に入力された式を D10 にコピーするか,D3:D10 の範囲ですべてのセルに D3 の式をオートフィ図 71 y として直線上の2 点を加えた様子ルしてもよい. 同じ直線上の点間を結んでも一直線に見えるはずである. 少なくとも D3 と D10 に傾き, 切片で定まる直線のデータがあるとして, そのデータ系列に名前を与えるため,D2 に文字列として y を入力した結果が図 71 である.

7 次に B2:D10 をアクティブ化し, グラ フ化して y の系列のみ線形近似した結 果が図 72 である. 表計算ソフトの基礎 37 図 72 y に関する 2 点間線形近似の描画結果 5.3 線形近似直線における最小二乗法図 73 の左に示した線形近似は, 最小二乗法によるものである. 一方, 図 73 の右に示した線形近似は, 他の考え方に基づく一例である. 直感的には, データ点と直線との距離が最も低くなる図 73 の右側の例の方が線形近似を行う際には, より良い選択肢であるように思える. 図 73 最小二乗法による線形近似 ( 左 ) と他の例 ( 右 ) 最小二乗法により, 図 73 左のようになるのは, 最小二乗法では, 線形近似直線とデータ間の距離を考慮する際, 図 74 のように, データ点から y 軸に沿って下ろした図 74 の赤い線 ( 誤差 )e 1 ~e 8 を距離とした距離総和を最小化するように傾きと切片を算出する手法を採用しているからである. 図 74 最小二乗法による線形近似曲線との差をとる手法

8 表計算ソフトの基礎 38 最小二乗法における線形近似は, 以下の z を最小化するような a と b を求める手法となっている. z= N {y k y k} 2 = k=1 k=1 N {y k (b+ax k )} 2 つまり, 最小二乗法とは, 実測値のうち,y k と推定値 y k の偏差の二乗和 z が最小になるような直線を求 める方法であると言える.z が最小値をとるのは z が極小値 ( 接線の傾き 0) をとるときなので,z を a と b でそれぞれ偏微分し, それらを 0 とおいて, 連立方程式を解いて a と b を求めることができる. ここでは途中を省くが, 連立方程式を解くと, 以下 のように,a と b の算出式を求めることができる. a の右辺に着目すると, N x k y k =x 1 y 1 +x 2 y 2 + +x N y N k=1 N N x k y k =(x 1 +x 2 + +x N ) (y 1 +y 2 + +y N ) k=1 k=1 N x 2 k=x 2 1+x x 2 N k=1 N 2=(x1 x k +x 2 + +x N ) 2 k=1 であり, 値は,x k, y k と N に基づいて表現され, 演算は, 四則演算と累乗で表現されていることから, データが得られると, 表計算の基本的な演算で容易に計算が可能である. 例えば, 新たなシートを用意し, データ D を,B 列にと C 列にコピーしたのち, 各 k ごとに,x k y k を D 列で,x 2 k を E 列で計算することにすると, 図 75 のようになる.

9 表計算ソフトの基礎 39 これは,D2 に = B2 * C2 E2 に = B2^2 を入力し,D2:E2 を同時にアクティブ化して,9 行目までオートフィルした結果である. また,N については,F2 に = ROWS(B2:B9) を入力し, これを演算に利用する. 図 75 最小二乗法で傾き a を表計算する例このような準備をすれば, 傾き a を計算することは, 容易である. 例えば, 空いているセルを使い, 分子と分母に分けて計算してみよう. 分子は, 例えばセル G2 に = $F$2 * SUM(D2:D9) - SUM(B2:B9) * SUM(C2:C9) 分母は, 例えばセル G3 に = $F$2 * SUM(E2:E9) - SUM(B2:B9)^2 と入力し セル G4 で, = G2 / G3 すると, セル G4 の計算結果は となる. これは,SLOPE 関数の結果と一致する. また, 最小二乗法による傾き a と切片 b の計算式を見比べると, 分母は同じであるから, 例えばセル G6 に切片 b の計算式について, セル G3 の分母の演算結果と, 分子の式を = (SUM(E2:E9) * SUM(C2:C9) - SUM(D2:D9) * SUM(B2:B9)) / G3 と入力すれば, となる. これは INTERCEPT 関数の結果と一致する. また,5.2.1 で, 線形近似直線をグラフへ挿入した際, 近似曲線の方程式を挿入して表示された結果 f(x) = x と SLOPE 関数や INTERCEPT 関数, さらには LINEST 関数の結果と同じ線形近似直線になっていることから,calc では, いずれの場合にも最小二乗法による計算法を採用していることがわかる. < 線形近似直線を求める際の最小二乗法とは異なるアプローチ> 図 73 右側のような線形近似を行う方法は, また別に存在する. ここでは詳しく紹介しないが, 図 76 のように線形近似曲線へ垂線を下ろす方法が一例である. 図 76 線形近似曲線に垂直な偏差に着目する方法

10 表計算ソフトの基礎 曲線近似 5.2.1に示したように, 表計算ソフトでは, データに対し, 近似曲線を挿入することができる機能があった. その中には, 線形近似直線だけでなく, 多項式近似, 指数近似, 累乗近似, 対数近似などいくつかの選択肢があった. ここでは, 参考として, 関数を用いて曲線をフィッティングさせる際の係数や定数項の値を算出する方法を例示する 多項式近似線形近似直線を求める際の SLOPE 関数による傾きや INTERCEPT 関数による切片の値と LINEST 関数が返す値は, 同じであった. ただし,LINEST 関数は 線形近似直線を求めるだけの関数ではなく, もともとはデータに基づき, 多項式近似を行うための関数である. つまり, f (x) = a x + b という直線近似 ( 一次式 ) だけでなく,m 個の係数 a 1,,a m と, 定数項 b による f (x) = a m x m + + a 2 x 2 + a 1 x + b のような多項式に基づいた曲線近似を行うこともできる. このとき, m 個の係数 a 1,,a m のうち,m 番目の係数を求めるには =INDEX(LINEST(y 値のセル範囲, x 値のセル範囲 ^{m, m-1,, 1}), 1, m ) のようにする. また, 定数項 b を求めるには, 以下のようにする. =INDEX(LINEST(y 値のセル範囲, x 値のセル範囲 ^{m, m-1,, 1}), 1, m+1) 指数近似 f (x) = a exp( b x ) = a e b x において, 係数 a は, =EXP( INTERCEPT( LN(y 値のセル範囲 ), x 値のセル範囲 ) ) または, =EXP( INDEX( LINEST(LN(y 値のセル範囲 ), x 値のセル範囲 ), 1, 2) ) として求めることができる. また, 定数項 b は, =SLOPE( LN(y 値のセル範囲 ), x 値のセル範囲 ) または, 次のようにして求めることができる. INDEX(LINEST( LN(y 値のセル範囲 ), x 値のセル範囲 ), 1, 1 ) 累乗近似 f (x) = a x b において, 係数 a は, =EXP( INTERCEPT( LN(y 値のセル範囲 ), LN(x 値のセル範囲 )) ) または, =EXP( INDEX( LINEST(LN(y 値のセル範囲 ), LN(x 値のセル範囲 )), 1, 2) ) として求めることができる. また, 定数項 b は, =SLOPE( LN(y 値のセル範囲 ), LN(x 値のセル範囲 )) または, 次のようにして求めることができる. INDEX(LINEST( LN(y 値のセル範囲 ), LN(x 値のセル範囲 )), 1, 1 )