見出し画像

数学とPython 7 マクローリン展開 ～可視化

ネイピア DS

はじめに

シリーズ数学とPythonのご案内

「シリーズ数学とPython」は、数学の学習中に「解けない、無理～」と焦ったときに、Pythonで数値の動きや可視化を行って、理解の糸口を見つけたときのことを記事にします。

データサイエンス数学ストラテジスト上級公式問題集

この記事は「データサイエンス数学ストラテジスト上級公式問題集」の問題の解読中に見つけたヒントを取り扱います。

今回取り組む問題

問題３７「関数のマクローリン展開－微分の達人になろう！」

Python実装

やりたいこと

問題で与えられた１変数関数$${f(x)}$$をマクローリン展開して、グラフを描画します。
実は、自分が解いた答案が、公式問題集に記載の解答・解説と合わなかったので、Pythonの力を借りて計算の検算を試みたのです。
その後、解答・解説の記載に誤植があったことを正誤表で見つけました（納得）。

作戦

sympy を用いて、公式問題集の出題内容（$${x^4}$$の係数）をサクッと確かめてから、グラフを描画します。

実装の開始

インポート

### インポート
import sympy
sympy.init_printing;

関数$${\boldsymbol{f(x)}}$$の定義

### 変数、関数の定義

# 変数xの定義
sympy.var('x')

# 出題の関数f(x)を定義
question_fx = sympy.sqrt(1 -x +2*x**2)
display(question_fx)
出力イメージ

関数$${\boldsymbol{f(x)}}$$をマクローリン展開する 
sympy の series （ベキ級数展開）でマクローリン展開を行います。

【関数$${f(x)}$$のマクローリン展開】
$${f(x)}$$の第$${n}$$次導関数を$${f^{(n)}(x)}$$とするとき、
$${f(x) \\ = f(0) + \cfrac{f^{(1)}(0)}{1!}x + \cfrac{f^{(2)}(0)}{2!}x^2 + \cfrac{f^{(3)}(0)}{3!}x^3 + \cdots + \cfrac{f^{(n)}(0)}{n!}x^n + \cdots \\ =\sum^{\infty}_{k=0} f^{(k)}(0) \cfrac{x^k}{k!} }$$

引数 n=5 と設定することで、$${x^4}$$まで展開します。
$${x^5}$$以降は剰余項$${O(x^5)}$$にまとまります。

### f(x)をマクローリン展開する　
answer_me = sympy.series(question_fx, x, n=5) # n：展開する指数+1
display(answer_me)
出力イメージ

$${x^4}$$の係数は$${-\cfrac{21}{128}}$$でした。

マクローリン展開の結果をプロット
sympy の plot でグラフを描画します。
removeOメソッドで余剰項を除去します。

### マクローリン展開の結果をプロットする removeO():剰余項O()を削除
sympy.plot(answer_me.removeO(), question_fx, legend=True, ylim=(-8, 8));
出力イメージ

ブルーの曲線がマクローリン展開で得た式です。
峰が２つあります。
魚が水面から顔を出している感じがします🦈

オレンジの曲線は元の関数$${f(x)}$$です。
ブルーの曲線の凹の部分（魚の口）に当てはまっています。

今回も可視化できてよかったです！
（数式と比べて優しい温もりを感じます😉）


参考：テーラー展開

sympy の series は 引数 x0 に$${a}$$の値を指定することで、テーラー展開を行えます。
出題の関数$${f(x)}$$を$${x=1}$$の周りでテーラー展開してみます。
引数 x0=1 とします。

### f(x)をテーラー展開する　f(x)=question_fx, x=a, a=x0=1
sympy.series(question_fx, x=x, x0=1, n=5)
出力イメージ


おまけ
Pythonで問題集の問題を解くコードです。
Jupyter Notebook のセル（計算ステップ）ごとに、計算結果を表示しています。

### インポート
import sympy
sympy.init_printing;
### 変数、関数の定義

# 変数xの定義
sympy.var('x')

# 関数f(x)の定義
f = sympy.sqrt(1 -x + 2*x**2)
f
### 微分を繰り返す
# １回目
f_d1 = sympy.diff(f, x)
f_d1
# ２回目
f_d2 = sympy.diff(f_d1, x)
f_d2
# ３回目
f_d3 = sympy.diff(f_d2, x)
f_d3
# ４回目
f_d4 = sympy.diff(f_d3, x)
f_d4
### ４階微分の結果を整理（因数分解）
f_d4_fact = sympy.factor(sympy.expand(f_d4))
f_d4_fact
### （おまけ）４階微分を一括して実行する方法 diffの第３引数でn階微分を指定
sympy.factor(sympy.diff(f, x, 4))
### 第4次導関数f_d4_factのxに0を代入
f_d4_fact_0 = f_d4_fact.subs(x, 0)
f_d4_fact_0
### f(^4)(0) / 4! を計算する　★解答
f_4_coef = f_d4_fact_0 / sympy.factorial(4)
f_4_coef

おわりに

自力の計算結果をsympyで検証することができました。

おわり

ブログの紹介

noteで３つのシリーズ記事を書いています。
ぜひ覗いていってくださいね！

１．のんびり統計
統計検定２級の問題集を手がかりにして、確率・統計をざっくり掘り下げるブログです。
雑談感覚で大丈夫です。ぜひ覗いていってくださいね。
統計検定２級公式問題集CBT対応版に対応しています。

２．Python機械学習プログラミング実践記
書籍「Python機械学習プログラミング PyTorch & scikit-learn編」を学んだときのさまざまな思いを記事にしました。
この書籍は、scikit-learnとPyTorchの教科書です。
よかったらぜひ、お試しくださいませ。

３．データサイエンスっぽいことを綴る
統計、データ分析、AI、機械学習、Pythonのコラムを不定期に綴っています。
「統計」「Python」「数学とPython」「R」のシリーズが生まれています。

この記事が気に入ったらサポートをしてみませんか？