91160


Pythonの全くの初心者で、分かりにくい質問かと思いますが、投稿させていただきます。  
時系列モデルのシミュレーションを行うためのコード作成に関する質問です。

具体的には、  
`X(t)=0.5*X(t-1)+e(t) for t<=50`  
`X(t)=X(t-1)+e(t) for t>50`  
（ただし、e(t)は、平均0、標準偏差1の正規分布の乱数）  
について、  
（１）： t=100のときのX(100)の値をシミュレートする関数のコードを作成する  
（２）：
適当な値X(0)を入力した場合のX(1),X(2),X(3),...X(100)の値をプロットしたグラフ（横軸にt=1,2,3,...100、縦軸にX(1),X(2),X(3),...X(100)）を作成する  
が目標です。

（１）については、以下のコードを作成しましたが、正しいでしょうか。

[code]

    import numpy as np
    x = 0
    for i in range(50):
        x = 0.5*x + np.random.randn()
        print(x)
    for i in range(50):
        x += np.random.randn()
        print(x)
    
[/code]

（２）については、例えば、numpyを利用する方法も試しましたが、関数（x = 0.5*x +
np.random.randn()）から導かれる計算結果を行列やベクトルに変換したり組み込む方法が分かりませんでした。

なお、  
Python 3  
Jupyter Notebook (anaconda 3)  
を使用しています。

かなり色々なサイトで探したのですが、類似のケースが見つからず、こちらに投稿させていただいた次第です。  
分かりにくい質問で申し訳ありませんが、何卒よろしくお願いいたします。


CC BY-SA 4.0
2022-09-18T07:47:50.153
91155
2022-09-19T06:00:22.580
2022-09-18T11:32:37.350
43025
54440
question
python
Pythonの時系列モデルのシミュレーション


`print(x)`で単に表示している所をグラフ作成用のリストに追加していく処理にすれば良いのでは？  
あと初期値X(0)を外から与えるパラメータとする関数として定義すれば良いでしょう。  
以下のように出来ると思われます。`####`でコメントした所が変更点です。

[code]

    import numpy as np
    
    #### X(0) の値をパラメータとするシミュレーション用Y軸データ作成関数を定義
    def simfunc(initial):
        x = initial          #### X(0) は外部から与えられたパラメータ値とする
        y = []               #### グラフのY軸データ用リスト初期化
        for i in range(50):
            x = 0.5*x + np.random.randn()
            y.append(x)      #### グラフのX軸1～50用 Y軸データリスト追加
        
        for i in range(50):
            x += np.random.randn()
            y.append(x)      #### グラフのX軸51～100用 Y軸データリスト追加
        
        return y
    
    import matplotlib.pyplot as plt #### 以下グラフ作成・表示
    
    y = simfunc(np.random.randn())  #### グラフのY軸データ作成関数呼び出し
    x = list(range(1,101))          #### グラフのX軸データリスト(1～100)
    
    plt.plot(x, y)
    plt.show()
    
[/code]


CC BY-SA 4.0
2022-09-18T13:28:31.003
91160
2022-09-18T15:31:19.617
2022-09-18T15:31:19.617
26370
26370
91155
answer


作図の都合上, 前半の数値を `0.99`で  
計算には以下を利用

  * (docs.python.org) [itertools.accumulate](https://docs.python.org/ja/3/library/itertools.html#itertools.accumulate)
  * NumPy [numpy.cumsum](https://numpy.org/doc/stable/reference/generated/numpy.cumsum.html)

[code]

    from itertools import accumulate
    import matplotlib.pyplot as plt
    import numpy as np
    
    fhalf = list(accumulate(np.random.randn(50), lambda a, b: a *.99 +b, initial=150))
    lhalf = np.cumsum(np.random.randn(50)) +fhalf[-1]
    
    fig = plt.figure()
    plt.plot(fhalf +lhalf.tolist())
    plt.savefig('plt.png')
    
[/code]

[![plt.plot
の結果](https://i.stack.imgur.com/gM7wf.png)](https://i.stack.imgur.com/gM7wf.png)

* * *

質問がなんらかの授業の課題のようにも見えたので, 回答付けるのどうかと判断迷ったけれど  
とりあえず付けてみました  
[https://ja.meta.stackoverflow.com/questions/1805/学校の宿題は回答するべきでしょうか](https://ja.meta.stackoverflow.com/questions/1805/%E5%AD%A6%E6%A0%A1%E3%81%AE%E5%AE%BF%E9%A1%8C%E3%81%AF%E5%9B%9E%E7%AD%94%E3%81%99%E3%82%8B%E3%81%B9%E3%81%8D%E3%81%A7%E3%81%97%E3%82%87%E3%81%86%E3%81%8B)


CC BY-SA 4.0
2022-09-19T04:05:20.810
91170
2022-09-19T06:00:22.580
2022-09-19T06:00:22.580
43025
43025
91155
answer