フーリエ変換で周期を発見 part2 | 科学アラカルト [実践] | 逆ポーランド電卓の実践ウェブ rpn hacks!

逆ポーランド電卓の実践ウェブ rpn hacks!

逆ポーランド電卓rpnの実践ウェブ   
rpn hacks! アールピーエヌ・ハックスサイトマップ

rpn | 実戦 | 科学アラカルト | 科学に関する話題をrpnで探求!調べて実験してみよう。

HOME > 実践 > 科学アラカルト > フーリエ変換で周期を発見 part2 hatena yahoo buzzurl livedoor del.icio.us nifty newsing twitter facebook rss ソーシャルブックマーク

フーリエ変換で周期を発見 part2

前のページに戻るLinkIcon

黒点の周期

 ここ数年、十数年は太陽の黒点数から目が離せない状態です。何世紀も前にあった小氷河期(マウンダリー極小期:1645年から70年間)の再来が迫っているかもしれないのです。

情報   太陽の黒点数に関しては、実践コーナー(科学アラカルト)に地球温暖化と太陽の黒点数があります。事前知識として読んでおくとよいでしょう。

今から遡ること300年、ガリレオの時代からの太陽観測のデータを使って、黒点数の推移を見ることができます。黒点の数は太陽活動の活発度を示す指標として知られており、黒点数が多ければ太陽は活発、少なければ太陽は低調です。

そして、黒点数の増減には10~13年の周期があることが知られています。その周期をデータから導き出してみましょう。果たして、フーリエ変換の計算で得られるのでしょうか。

ちなみに、この黒点数の周期は1843年には知られていたそうです。フーリエ没後13年ですね。

まず、1700年から300年間にわたる年間平均の黒点数データが"sunspot.txt"にあるとします。このデータを時系列に並べてグラフにすると次のようになります。

  >rpn 1 -c rownum <sunspot.txt | xyp -x,300 -y,200 -s50,100 -m -w60
  ^y 200                                            *
  |                                                 *
  |
  |              *                                *  *   * *
  |                          *     *              ***    * *
  |    *         * *        *  *
  |      *         *               *            **   *   *
  -    * *     *             * * *  *       *         **  **
  |      * ** * ** *         *          *    *   ***    *    *
  |    * **    *  *        *     * *   **     **** *   *
  |* *     * ** *   *      ******** *** * * ****     * * * **
  |  *  *   ***     * *  * *  *   * * * * *** *  *    *   *
  ** ************ *   ***** **  * ** * * *   * ******* ** * *
  **  * * **  *  ** **  ***   * *** * * * ** *** * *  * * * *x
  *** * * **** ** *  * **** * * *  ***** * ** *** * * * *  3*0
  +-*------|---------|-*-------|---------|---------|--------->


プロットが密集しているのでよく分かりませんが、黒点数の最大は200程度であることが分かります。最初の50年間だけを拡大表示してみましょう。

  >rpn 1 -c rownum <sunspot.txt | xyp -x,50 -y,200 -s50,100 -m
  ^y 200
  |
  |
  |
  |
  |                    *
  |                             *
  -                     *        *
  |                            *         *
  |                    * *    *  *
  |   *         *                       *
  |            *          *
  |  * *      *  **   *   *   *   *    *
  | **  *          * *       *     ** *  x
  **     ** **      *      **        *  50
  +--------**---------------------------->


確かに周期がありそうですね。50年間に5回のピークなので約10年の周期を持っていそうです。

黒点数の増減周期をフーリエ変換で見つける

 では、このデータを使ってフーリエ変換後、スペクトルを出してみましょう。フーリエ変換はrpnでは1行でOKでしたね。

  >rpn -c fourier <sunspot.txt >specspot.txt


specspot.txtにスペクトルのデータが格納されているはずです。グラフ化してみましょう。

  >rpn 0 -c rownum <specspot.txt | xyp -x,299 -y,5000 -s50,1000 -m
  ^y 5000
  |
  -  *                               *
  |
  |  *                               *
  -
  |
  *                                     *
  -  *                               *
  |  *                               *
  *  **                             **  *
  *   *                             *   **
  |*******                       ********
  *********                     *********x
  |**  ****************************** **99
  +-----|------|------|-----|------|----->


綺麗なスペクトルが出てきました。例のごとく150の点がナイキスト周波数で、後半は無視します。加えて、最初の数値を除いた前半部での最大値が周期の候補になります。前半部分を拡大してみましょう。

  >rpn 0 -c rownum <specspot.txt | xyp -x,49 -y,5000 -s5,1000 -n -m
  ^y 5000
  |
  4000                 *
  |
  |                      *
  3000
  |
  | *
  2000               *
  |                     *
  |   *                  *     *
  **0* *                     *
  |     **   *    * *     *  **  *
  |  *   ***    ** * **    **   * * **   *
  |o        ****                 * * *****
  +---5---10--15--20--25--30--35--40--45->  ^y 5000


横軸が25から30までの間にピークが見られます。実際にスペクトルの値が大きい順に並べて、数値で確認してみましょう。

  >rpn 0 -c rownum <specspot.txt | rpn 1 149 -c lookup | rpn x | sort -n
     :
   (上略)
     :
  1930.98 28
  2316.63 25
  2576.26 3
  3635.81 30
  4000.54 27


最大は27ですね。300年間で27番目なので、1周期が何年に相当するか計算します。30番目と3番目も一緒に。

  >rpn 300 27 / 300 30 / 300 3 /
  11.1111 10 100


結果は11.1年でした。上記の説明とほぼ合ってますね。次の候補は10年周期ですし、上位TOP5に25、27、28、30とあるので、10~12年くらいの周期を示しています。

最後にですが、上位の3番目にもピークが出ています。計算すると100年周期です。もしかすると黒点数の周期には11年前後の周期と100年前後の周期が混じっているのかもしれません。

コレログラムでの黒点数検証

 ちなみに前述の記事でも使用したコレログラムという指標を使って黒点数の推移を分析すると次のようになります。

  >rpn -c crlogram <sunspot.txt | rpn 0 -c rownum | xyp -x,300 -m
  *y 1
  *
  |
  |*
  ***
  | * *
  |*  **      **** *     ****** *
  * ** ***********************************
  +*-*************************************
  *o*****       **** *   *  ********     x
  *****
  *
  |
  |
  |
  |-1


前半に相関係数が高く、後半に向かって急速に値が減っていきます。ただ、興味深いことに100年、200年くらいの部分でなだらかなピークが見えるような気もします。まず、前半の50年分だけを拡大してみましょう。何かの相関が見られるでしょうか。

  >rpn -c crlogram <sunspot.txt | rpn 0 -c rownum >tmp
  >xyp -x,50 -s5 -n -m <tmp
  *y 1
  *
  |
  |      **
  |*     * *      **
  |              * *      **
  |     *   *             * *    ****
  | *           *   *    *           *  5*
  +--5---10-*15--20--25-*30--*5--*0--*5--*
  |o   *        *    *        * *     ***x
  |  **      ***      **      **
  |  *
  |
  |
  |
  |-1


見事に自己相関の傾向が見えます。1つ目の山は無視するとして、10年、11年あたりにピークが出ています。そのn倍年に同様にピークが出ていることから、最小の周期が10、11年ということが分かりますね。

ただし、上記のグラフではフーリエ解析で出てきた100年の周期は判別はできません。50年から150年までの部分を拡大すると何か見えてくるでしょうか。

  >xyp -x50,150 -y-.5,.5 -s10 -n -m -w60 <tmp
  ^y 0.5
  |
  |
  |
  |                           *     **     *
  |                           **    **    *      *
  |**                  ***   *  *     *     *    **     *
  *0*    ***    ***    *  *  *     *      *     * *    * * 150
  +--*-6**--*70**--**-*--9*----1*0-*-11*---12*---13*---*40*--*
  *   **     **     **     **    **      *   *  *     *   ****
  *                                    **     **    **     *
  |                                           *      *
  |
  |
  |
  |-0.5


僅かにですが、100~120年のあたりに相関係数の盛り上がりが見えます。感覚的には偶然とは思えないような気はします(2倍の200年の部分にも盛り上がりがあるので)。

ただし、他人への説明(説得)は難しいでしょう。その点、フーリエ解析の結果なら説明しやすいかもしれません。このようにコレログラムでもフーリエ解析の結果を支持できることが分かりましたね。

株価の周期性の探索

 株価には周期があるのかないのか。永遠の疑問です。周期といっても短いもの(数日単位)から長いもの(数年・数十年)まであります。経済政策、税制改正、為替変動、投機筋の資金流入、投資機関によるリバランス…。

仮に株価に周期があるとしても、その要因は数え上げればキリがありませんし、分析もほぼ不可能です。でもフーリエ解析なら何かの周期を炙り出せるかもしれません。

さっそく、日経平均株価を使って検証してみましょう。ただ、日毎の株価では細かすぎるので、月末の株価を使うことにします。期間はバブルや暴落を含む1985年から2010年までの26年とし、合計で312点のデータを集めます。なお、日経平均のデータはファイルの"nikkei.txt"にあるとします。

  >rpn 1 -c rownum <nikkei.txt | xyp -x,312 -y,40000 -s12 -m -w60 | npd
  ^y 40000   * 日経最高値
  |         **
  |        ****
  |       ** **
  |      **
  |    ***    ****
  |    **      ****        **   ITバブル       ライブドアショック
  |  **        *  * ***** *****     *                 リーマンショック
  |  **           *** * ***  **** ****           ***** (サブプライムローン)
  | **             *     *     ****  ***         **   **
  *** ブラックマンデー                      ****   ****     **   *
  |プラザ合意                            *****           ******
  |                                        *            *
  |   バブル景気     失われた10年            イザナミ景気        x
  |o                                                       312
  +-|-|--|-|-|--|-|--|-|-|--|-|--|-|-|--|-|--|-|-|--|-|--|-|->
  1985      1990        1995        2000        2005      2010


株価の推移をグラフにすると次のようになります。

この期間はプラザ合意での円高ドル安誘導に始まり、好調な経済(不動産投資全盛)を背景に株価が4万円に届こうかというバブル絶頂を経て、失われた10年を経験。その後、2000年前後のITバブルとその崩壊、BIS規制よる国内銀行の再編、郵政事業の民営化、米国サブプライムローン問題をきっかけとしたリーマンショック…とイベントが目白押しの四半世紀でした。

やはり、株価のグラフからは一見、これといった周期は感じられません。では、黒点数の分析で効果を発揮したコレログラムをまず試してみましょう。

  >rpn -c crlogram <nikkei.txt | rpn 0 -c rownum | xyp -x,312 -m
  *y 1
  **
  |*
  | *
  | **
  |  **
  |   ** *****
  |     **   **                      ****2
  +-----------*****-------****----***---->
  |o              ***    **  ******      x
  |                 *****
  |
  |
  |
  |
  |-1


コレログラムでも、これといったピークはないようです。直近の株価が影響するということ以外には目新しいものはありません。精々、150ヶ月頃(12年)に僅かに逆相関が見える程度です。

株価の増減周期をフーリエ変換で見つける

 それではフーリエ解析はどうでしょうか。目に見えない周期を炙り出してくれるでしょうか。まずはフーリエ変換とスペクトルの計算です。

  >rpn -c fourier <nikkei.txt >specnik.txt


次にスペクトルをグラフ化するために、最初の数値以外の最大値を求めます。

  >rpn 0 -c rownum <specnik.txt | rpn 1 311 -c lookup | rpn _, -c max
  1.04214e+06


最大は1.04e+6であることが分かります。縦軸の最大を1.5e+6にしてグラフにしてみます。

  >xyp -x,311 -y,1.5e+6 -m <tmp
  ^y 1500000
  |
  |
  |
  *                                      *
  |
  |
  *                                     *
  |
  |
  *                                     *
  |
  *                                     *
  **                                   **x
  ***************************************1
  +-------------------------------------->


典型的なグラフ形状のように見えますが、最初と最後のピークの立ち上がりが強すぎるので片対数で再度、表示してみます。

  >xyp -x,311 -yl,1.5e+6 -m <tmp | npd
  *                                     **
  *                                     *
  **                                   **
  ******                           ******
  |*********** ****     ****  **********
  |* ***********************************
  |      * ***** *************** *
  |    **     * **       ** *     **
  ~
  |
  +-------------------------------------->


やはり、良く見るスペクトルの形です。では、スペクトルの最初の方の値を見てみましょう。

  >rpn 10 -c head <tmp
  1 1.04214e+06
  2 422921
  3 702030
  4 281509
  5 191420
  6 83441.1
  7 155217
  8 155870
  9 6552.38
  10 120361


1番目の値が1桁違いで大きいことが分かります。あとは急激に小さくなっています。どうも、20番目くらいまでを拡大表示すれば十分そうです。

  >xyp -x,20 -y,1.5e+6 -s1 -n -m -w60 <tmp
  ^y 1500000
  |
  |
  |
  | *
  |
  |
  |       *
  |
  |
  |    *
  |
  |          *
  |             *     *  *     *        *                    x
  |o               *        *     *  *     *  *  *  *  *  * 2*
  +-1--2--3--4--5--6--7--8--9--10-11-12-13-14-15-16-17-18-19->


これといったピークは感じられませんが、あるとすれば1番目と3番目くらいでしょうか。

  >rpn 312 1 / 312 3 /
  312 104


312ヶ月は26年、104ヶ月は8.7年に相当します。ただ、26年周期に意味は感じられず(元々が26年間のデータ)、8.7年の周期(26年間に3回の繰り返し)にも理由を見つけ出すことは困難そうです。

四半世紀の日経平均を取り込み、その期間にあったバブルも暴落も全て飲み込んで、何の周期も見出せないことになります。そもそも株価はランダムウォークに例えられるくらい動きが不規則で、秒や分、時間、日、週、月、年のどの時間軸にもランダム性が色濃く出てきます。

今回の解析で月単位での周期性は見られませんでしたが、時間軸の違いによっては周期性が発見できるのかもしれませんね。

実践科学アラカルトに戻るLinkIcon

警告rpnプログラムを実行するには、rpn試用版かrpn標準版が必要です(バージョンの違いはこちら)。

情報pasteは講座サポートで公開されています。rownumはカンタン分析パッケージに同梱されています。lookupはユーティリティーパッケージに同梱されています。statinfo, headはrpnマイスターパッケージに同梱されています。crlogramはビジネス統計(トレンド編)に同梱されています。xypとnpdはrpnの姉妹ソフトウェアです。詳しくはプロダクトを参照ください。

科学アラカルト

調べて実験してみよう

※実践コーナーのTOP

紹介 rpnの利用シーンはこちら…

講座初めての人のrpn基礎もどうぞ
講座しっかり学べるrpn入門もどうぞ
講座すぐに使えるdos入門もどうぞ
物語データを見抜くojt物語もどうぞ

実践他の分野への挑戦は実践TOP

応用rpnアプリケーションは応用TOP

part2

調べて実験してみよう

※実践コーナーのTOP

講座初めての人のrpn基礎もどうぞ
講座しっかり学べるrpn入門もどうぞ
講座すぐに使えるdos入門もどうぞ
物語データを見抜くojt物語もどうぞ

実践他の分野への挑戦は実践TOP

応用rpnアプリケーションは応用TOP

書籍紹介

記事に関連した書籍

本ウェブサイトで扱った話題に関連した書物で、スタッフが実際に読了したものを紹介。

書籍数学の書籍
数の世界は思ったよりもエキサイティング。

  • 書籍統計の書籍
  • ビジネスで統計が使えるととっても有利。

書籍投資の書籍
失敗しない投資には広範囲で実践的な知識が必要。

警告バックスラッシュはエンマークに読み替えてください( IEのみ)。
バックスラッシュとエンマーク

警告文字で作られた図表や式が崩れることがあります。ブラウザによっては固定幅フォントをMSゴシックにするときれいに表示されます。それでも崩れる場合は図表や式をメモ帳にコピー後、閲覧下さい。

警告rpn試用版と標準版(2kリビジョン)はダブルクォートで囲ってください。

rpn 1 2 + ⇒ rpn "1 2 +"
rpn 1 -c foo ⇒ rpn "1" -c "foo"

ダブルクォートは省略できることが多いのですが、慣れない間は囲んだほうが無難です。なお、本ウェブサイトの記事ではrpn標準版(98リビジョン)を使用しているため囲っていません。詳しくは技術サポートの「rpn TIPS参照ください。

注意rpnの障害情報はこちら

警告rpn試用版の場合、複雑なプログラムや処理時間のかかるプログラムの一部には動作しないものがあるかもしれません。あくまで無料提供であることを勘案・了承ください。rpn標準版は、すべてのプログラムが動作します。