#author("2021-04-09T20:24:06+09:00","default:LabMember","LabMember") #author("2021-04-09T20:25:18+09:00","default:LabMember","LabMember") [[PythonTutorial]] * Pythonを使ったデータ解析その1 [#p8db5abe] #contents ** 準備 [#d44aa4e0] ex1.zipをダウンロードして展開し,ex1.ipynbを開く [[ex1.zip>https://granite.phys.s.u-tokyo.ac.jp/wiki/Lab/?plugin=attach&pcmd=open&file=ex1.zip&refer=PythonTutorialExample1]] をダウンロードして展開し,ex1.ipynbを開く *** 中身について [#y6c67512] - ex1.ipynb: メインのデータ解析をするnotebook - data1_*: 伝達関数の測定データ. 周波数/ゲイン[dB]/位相[deg]/コヒーレンス の順に列が並んでいる ** データの読み込み [#le086ee2] まず解析するデータの読み込みを行う.まずは numpy.loadtxt を使う. メインのデータはどのファイルも同じだが,拡張子やエンコードなどが微妙に違い,うまく対応したコードでないと読み込めない. これらのデータ形式は実際に実験室で測定した場合に得られる可能性のあるものを選んでいる. *** data1_1.txt [#la911823] ただのテキストファイル. np.loadtxt(読み込むファイルのパス) で読み込み, np.ndarray 形式で出力する. *** data1_2.txt [#r9173e0a] ↑のデータにヘッダーがついていて,カンマで区切られている. skiprowsでヘッダーの分だけスキップする行数を指定し,delimiter で区切り文字を "," に指定する. np.loadtxt(読み込むファイルのパス, skiprows=読み飛ばす行数, delimiter=区切り文字) *** data1_3.csv [#ne7b4ea0] 中身はdata1_2.txt と同じだが,ファイル形式がcsvになっていて文字コードがshift-jisになっている(普通はutf-8). ファイル形式が.csvだろうと.txtだろうと numpy.loadtxt で読み込めるが,問題は文字コード. np.loadtxt(読み込むファイルのパス, skiprows=読み飛ばす行数, delimiter=区切り文字, encoding=文字コード) 文字コードでエラーが出た時は大体shift-jisになっているせいだと思えばよい. *** data1_4.csv [#l026e577] これまで numpy.loadtxt で読んできたが,ここでは pandas.read_csv を用いる. 今回はデータ量が少ないために気にならないが,実は numpy.loadtxt の読み込みはめちゃくちゃ遅い. 対して pandas.read_csv はかなり高速(一説には100倍以上の速度)で読み込める. pd.read_csv(読み込むファイルのパス, skiprows=読み飛ばす行数, delimiter=区切り文字) 出力結果は pandas.dataframe という形式でそのままでは使いずらい. なので,numpy.ndarray 形式に変換して使うとよい. そのためには,pandas.dataframe の後ろに .values をつける. dataframe = pd.read_csv(読み込むファイルのパス, skiprows=読み飛ばす行数, delimiter=区切り文字) data = dataframe.values ** プロット [#tc4fb608] 読み込んだデータのプロットをする. この辺は個人の好みが出るところなので,自由にカスタマイズしてよし. *** プロットサイズの指定 [#d2c3357f] デフォルトでプロットすると,かなり小さな図になってしまう.プロットサイズを大きくするには plt.figure(figsize=プロットサイズ) と最初に指定する. デフォルトは[6.4, 4.8]だが,例えば[12,10]とするとほぼ倍のサイズになる. plt.figure(figsize=[12,10]) *** 対数スケール [#z81de0df] デフォルトでは線形スケールだが,両対数・片対数でプロットすることは多い. そういった場合は plt.xscale("log") # x軸を対数スケールに plt.yscale("linear") # y軸を線形スケールに などとすると,各軸のスケールを線形/対数に設定できる. またプロットの時点で plt.loglog(x,y) # 両対数プロット plt.semilogx(x,y) # x軸が対数,y軸が線形の片対数プロット とすることでも可能. *** 軸設定 [#c495aa3d] 軸ラベルは plt.xlabel(ラベル名) # x軸のラベル plt.ylabel(ラベル名) # y軸のラベル で設定できる. 軸の目盛はデフォルトでは自動で設定されるが,自分で指定する時は, plt.yticks(目盛のリスト) とする.~ 位相をプロットする時は90º刻みなどで指定するとわかりやすい. *** グリッド [#u67fed4c] グリッドは plt.grid() で描くことができる. デフォルトでは主軸のグリッドだけ引かれるが, plt.grid(which="minor") とすると副軸のグリッドを引ける. *** 複数のプロットを並べる [#k55de647] 複数のグラフを同じ画面上に配置したい時は, plt.subplot を使う. ax1 = plt.subplot(行数,列数,何番目か) プロットする時は ax1.plot(x,y) などとすればよい. 各プロットの幅を変えたい場合, matplotlib.gridspec を用いる. GS = matplotlib.gridspec.GridSpec(行数,列数) 例えば2つのグラフを縦に並べて,高さを 2:1 にしたければ GS = matplotlib.gridspec.GridSpec(3,1) ax1 = plt.subplot(GS[:2]) ax2 = plt.subplot(GS[2]) のようにすればよい. ** フィッティング [#ycc6f416] フィッティングをする方法はいくつかあるが,ここでは scipy.optimize.curve_fit を使う.~ 以下では,データに対して未知パラメーターが2つ(a,b)の関数f(x,a,b)をフィットさせたいとする. *** フィッティングする関数の定義 [#y7b94f57] フィットさせたい関数を定義する. def f(x, a, b): # 定義文 ここで引数の1番目は必ずデータのx軸の値になるようにする. 2番目以降の引数がフィットさせるパラメーターに対応する. *** フィッティング [#u1d57f34] scipy.optimize.curve_fit は scipy.optimize.curve_fit(関数名, データx, データy) と使う.出力結果は[推定した値,推定値の共分散行列]となる.~ 例えば,データ列 x, yに対して関数f(x,a,b)をフィッティングしてa, bを推定する場合には params, cov = scipy.optimize.curve_fit(f, x, y) などと書く. ここでは params[0] がaの推定値, params[1]がbの推定値で, sqrt(cov[0,0])がaのフィッティングの誤差, sqrt(cov[1,1])がbのフィッティングの誤差になる. *** 初期値を指定してフィッティング [#v15edeae] フィッティングの精度が初期パラメーターに敏感な場合はよくある.そういった時に,予め当たりをつけてからフィッティングするとうまく求まる. 初期パラメーターの指定は scipy.optimize.curve_fit(関数名, データx, データy, p0=(初期パラメーターのリスト)) とすればよい.