えーウイルスも落ち着きの傾向を見せていますが、このままだと嬉しいことに自粛解除になってしまいそうです。
この腐りきった肉体と精神のまま、自粛解除されても困ってしまいますよ。
適度な運動、適度なコミュニケーション、適度な仕事、私が求めているのはそういう生活です。笑
人混み嫌です。満員電車嫌です。疲れ切った人間関係嫌です。人類補完願います。
さてさて、こういう機会じゃないと、ブログの更新も滞ってしまうので、ネタ作りのためにもnumpyで遊んでみます。
まぁ、そうは言ったものの、やる気も起きず、FFTを軽くやります。”軽く”ね。
数学的なことは、やりません。面倒なのでね。
したがって、マジメに勉強したい方は、教科書でも買ってください。いくらでもあると思います。
信号処理はどんな分野でも必要な知識なので、勉強しておいて損はしないと思います。（求められる程度は分野によって違うでしょうけど）
では、課題設定から。
「よくわかんねーけど、ウイルスの新規感染者数のデータをフーリエ変換する」
です。
大事なポイントはよくわかんねーけどですね。
このラフさが僕の長所です。
昨今のウイルスの日ごとの感染者の数が出ているので、そちらをこちらからダウンロードします。
残念ながら、直接は使い物にならないので、データの整形をします。（感染者一人に対してのデータなので、日にちごとの総数がわからないため。）
ちょこっとなので、エクセル内でやってしまいます。
新規シートに2020/01/24からダウンロードした日の一日前までの日付をA列に作成します。一日作ったらバーッとコピーすれば作れます。
B列一番上には
=countif(シート１のE列,A1)
と入力して、以下コピーします。（シート１のE列はちゃんと選択してくださいね笑）
これでできた日ごとの新規感染者数のデータをコピーして新規でCSVで保存。
適当にエクセルでグラフを作れば、
になりますね。ワイドショーで見たやつだ！！！ってなると思います。
それにしても減ってますね。小並感
じゃあ、あとはFFTのプログラムを張り付けて完成。
import numpy as np
import matplotlib.pyplot as plt
# サンプリング周波数（day）
rate = 1
# データの読み込み
signal = np.loadtxt("kansensya_day.csv",delimiter=",",usecols=1)
#時間軸のグラフ
plt.plot(signal)
plt.xlim(0, signal.size)
plt.xlabel('time [day]', fontsize=20)
plt.title('Time series of signal', fontsize=20)
plt.show()
# パワースペクトル
p = np.abs(np.fft.rfft(signal))
# 各成分の周波数を返す
f = np.fft.rfftfreq(signal.size, d=1./rate)
#周波数軸のグラフ
plt.xlabel("Frequency [1/day]", fontsize=20)
plt.ylabel("Power spectrum", fontsize=20)
plt.yscale("log")
plt.plot(f, p)
plt.show()
サンプリング周波数は1日置きのデータなので、1/dayになります。
データは先程作った同ディレクトリのkansensya_day.csvから引っ張ってきます。
日付の情報はいらないので、usecolsで二列目を選択します。（0から数えますよ）
あとはフーリエ変換して、グラフにしてあげます。
まぁ、細かいところに注意しながらフーリエ変換するものですが、注意しないでフーリエ変換していきます。
数学的な厳密さは考えません。
できたグラフがこちら。
ふむふむ、なるほど、さっぱりだ。
まあ、こんな結果だと思ってましたが、見たかった一週間の周期がうっすら見えてるので、いいんじゃないでしょうか。
え、一周間の周期はどうやって見ればいいかって？
横軸は1/dayなので、1/7（=0.14くらい）のところを見ればいいんですよ。
縦軸は相対値だと思って、数字に意味はないと思った方がいいですね。
結論から言えば、fftをしてあげる相手が間違っているのですが、fftの練習にはちょうど良かったと思います。
全国各地のデータでやってみるのも面白いかもしれませんね。やりませんけど。