以下に正規分布の散布図と楕円(3σ)，最小二乗法のプロットを示します．楕円を描くメソッド confidence_ellipse() は，Plot a confidence ellipse of a two-dimensional dataset 記載されているものを使いました．
# -*- coding:utf-8 -*-
import numpy as np
from numpy.random import multivariate_normal as normal
import matplotlib.pyplot as plt
from matplotlib.patches import Ellipse
import matplotlib.transforms as transforms
# -------------------------------------------------------------------------
# https://matplotlib.org/devdocs/gallery/statistics/confidence_ellipse.html
# 入力引数
# x, y: 入力データ 配列, shape (n, )
# ax: matplotlib.axes.Axes
# n_std: float 標準偏差の n_std 倍のプロットを作成
# **kwargs: `~matplotlib.patches.Ellipse` が引き継ぐ
# 戻り値
# matplotlib.patches.Ellipse
# -------------------------------------------------------------------------
def confidence_ellipse(x, y, ax, n_std=3.0, facecolor='none', **kwargs):
''' 分布の楕円を作成する． '''
if x.size != y.size:
raise ValueError("x and y must be the same size")
cov = np.cov(x, y)
mean_x, mean_y = np.mean(x), np.mean(y)
pearson = cov[0, 1]/np.sqrt(cov[0, 0] * cov[1, 1])
ell_radius_x = np.sqrt(1 + pearson)
ell_radius_y = np.sqrt(1 - pearson)
ellipse = Ellipse((0.0, 0.0),
width=ell_radius_x*2, height=ell_radius_y*2,
facecolor=facecolor, **kwargs)
scale_x = np.sqrt(cov[0, 0]) * n_std
scale_y = np.sqrt(cov[1, 1]) * n_std
transf = transforms.Affine2D() \
.rotate_deg(45) \
.scale(scale_x, scale_y) \
.translate(mean_x, mean_y)
ellipse.set_transform(transf + ax.transData)
return ax.add_patch(ellipse)
if __name__ == "__main__":
# --- プロットデータの作成 ---
xlim, ylim = [-10, 10], [-10, 10]
mean = [-1, -2]
cov = [[8, -3], [-3, 4]]
x, y = normal(mean, cov, 1000).T
conf = np.polyfit(x,y,1)
xls = np.linspace(xlim[0], xlim[1], 256)
yls = conf[1]+conf[0]*xls
# --- プロットの作成 OO-API ---
fig = plt.figure()
ax = fig.add_subplot(1,1,1)
ax.set_title("2D normal distribution", fontsize=20, fontname='serif')
ax.set_xlabel("x data", fontsize=20, fontname='serif')
ax.set_ylabel("y data", fontsize=20, fontname='serif')
ax.tick_params(axis='both', length=10, which='major', labelsize=16)
ax.tick_params(axis='both', length=5, which='minor')
ax.set_xlim(xlim)
ax.set_ylim(ylim)
ax.minorticks_on()
ax.grid(b=True, which='major', axis='both')
ax.scatter(x, y, s=5, marker='o', edgecolor="b", label='scatter data')
confidence_ellipse(x, y, ax, n_std=3,
label=r'$3\sigma$', edgecolor='red')
ax.plot(xls, yls, 'g-', label='Least squares')
# ----- プロットの表示 -----
ax.legend()
plt.show()
# ----- プロットファイルの作成 -----
fig.savefig('normal_dist.pdf', orientation='portrait')
fig.savefig('normal_dist.svg', orientation='portrait')
fig.clf()