Go to Qiita Advent Calendar 2024 Top
LoginSignup
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 3 years have passed since last update.

@ochimusha01(yasunori matuki)

【数理考古学】誤差関数(ERF=Error Function)と相補誤差関数 (ERFC=Complementary Error Function)。

Last updated at Posted at 2020-11-13

【初心者向け】「三角不等式(Triangle Inequality)」の体感方法？ - Qiita
image.gif
image.gif
image.gifimage.gif
正規分布(Normal Distribution)概念の根底にある中心極限定理（CLT=Central Limit Theorem）をガウスはこう解釈しました。
【無限遠点を巡る数理】無限遠点としての正規分布と分散概念の歴史

  • 小さな誤差ほど比較的観測され易い。
  • 大きな誤差ほど比較的観測され難い。
  • 従って、一定以上の誤差を切り捨てる事によって観測精度(Observation Accuracy)を確定する事が出来る。

その具体的表現(Expression=数式化)が誤差関数(ERF=Error Function)と相補誤差関数 (ERFC=Complementary Error Function)となります。
誤差関数 - Wikipedia
Rで誤差関数・相補誤差関数erf,erfcを使う

#誤差関数(ERF=Error Function)
誤差関数(ERF=Error Function)erf(x)=2/sqrt(π)×integrate(exp(-t^2), 0, x) 
image.png

erf <- function(x) 2 * pnorm(x * sqrt(2)) - 1
plot(erf,xlim=c(-3,3),ylim=c(-1,1),main="Error Function(ELF)",xlab="x",ylab="ELF(x)",col=rgb(1,0,0))
abline(h=0)
abline(v=0)

x軸とy軸を逆転させた逆関数(Inverse Function)はこんな感じになります。
image.png

erf <- function(x) 2 * pnorm(x * sqrt(2)) - 1
cy<-seq(-3,3,length=61)
cx<-erf(cy)
plot(cx,cy,type="l", xlim=c(-1,1),ylim=c(-3,3),main="Inverse ELF",xlab="x",ylab="Inverse ELF(x)",col=rgb(1,0,0))
abline(h=0)
abline(v=0)

これを以下の投稿で使用した「円形尺(Circular Scale=数直線概念の延長線上に現れる角度概念導入以前の曲座標系)」に射影(Projection)するとこんな感じになります。
【初心者向け】「観測者問題」のあっけない解決方法? - Qiita
image.png


erf <- function(x) 2 * pnorm(x * sqrt(2)) - 1

erf_scale00<-function(inbetween){

c0<-seq(0,pi*2,length=61)
cx<-cos(c0)
cy<-sin(c0)

c0<-seq(0,3,length=31)

for (i in c0){
plot(cx*(erf(i+inbetween)),cy*(erf(i+inbetween)),type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="",xlab="",ylab="",col=c(200,200,200,200))
par(new=T)
}

plot(cx,cy,type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="Error Function(Elf) Scale",xlab="Cos(θ)",ylab="Sin(θ)i", col=rgb(0,1,0))

abline(h=0,col=rgb(1,0,0))
abline(v=0,col=rgb(1,0,0))

# 凡例を書き添える 。
legend("topright", legend=c("Circle(Radius=1,Augment=θ)","x=y=0"), lty =c(1,1),col=c(rgb(0,1,0),rgb(1,0,0)))
}

erf_scale00(0)

アニメーション化(0→∞、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/30,2/30,0,1/30,2/30,0,1/30,2/30,0,1/30,2/30)
saveGIF({
for (i in Time_Code){
　erf_scale00(i)
}
}, interval = 0.1, movie.name = "erf_scale_back.gif")

アニメーション化(∞→0、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,2/30,1/30,0,2/30,1/30,0,2/30,1/30,0,2/30,1/30)
saveGIF({
for (i in Time_Code){
　erf_scale00(i)
}
}, interval = 0.1, movie.name = "erf_scale_forward.gif")

アニメーション化(どちらにも見える、中割り一枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/20,0,1/20,0,1/20,0,1/20,0,1/20,0,1/20)
saveGIF({
for (i in Time_Code){
　erf_scale00(i)
}
}, interval = 0.1, movie.name = "erf_scale_even.gif")

そのまま水平線/地平線(Horizon)を表すSin(θ)を採用するとピッチが倍に見えます。
Rplot17.png

sin_scale <- function(x) Im((0+1i)^(2*x))

sin_scale00<-function(inbetween){

c0<-seq(0,pi*2,length=61)
cx<-cos(c0)
cy<-sin(c0)

c0<-seq(0,1/2,length=31)

for (i in c0){
plot(cx*(sin_scale(i+inbetween)),cy*(sin_scale(i+inbetween)),type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="",xlab="",ylab="",col=c(200,200,200,200))
par(new=T)
}

plot(cx,cy,type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="Sine Scale",xlab="Cos(θ)",ylab="Sin(θ)i", col=rgb(0,1,0))

abline(h=0,col=rgb(1,0,0))
abline(v=0,col=rgb(1,0,0))

# 凡例を書き添える 。
legend("topright", legend=c("Circle(Radius=1,Augment=θ)","x=y=0"), lty =c(1,1),col=c(rgb(0,1,0),rgb(1,0,0)))
}

sin_scale00(0)

アニメーション化(0→∞、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/180,2/180,0,1/180,2/180,0,1/180,2/180,0,1/180,2/180)
saveGIF({
for (i in Time_Code){
　sin_scale00(i)
}
}, interval = 0.1, movie.name = "sin_scale_back.gif")

アニメーション化(∞→0、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,2/180,1/180,0,2/180,1/180,0,2/180,1/180,0,2/180,1/180)
saveGIF({
for (i in Time_Code){
　sin_scale00(i)
}
}, interval = 0.1, movie.name = "sin_scale_forward.gif")

アニメーション化(どちらにも見える、中割り一枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/120,0,1/120,0,1/120,0,1/120,0,1/120,0,1/120)
saveGIF({
for (i in Time_Code){
　sin_scale00(i)
}
}, interval = 0.1, movie.name = "sin_scale_even.gif")

これを是正するにはピッチを半分にするのです。
Rplot18.png

sin_scale <- function(x) Im((0+1i)^(2*x))

sin_scale01<-function(inbetween){

c0<-seq(0,pi*2,length=61)
cx<-cos(c0)
cy<-sin(c0)

c0<-seq(0,1/2,length=16)

for (i in c0){
plot(cx*(sin_scale(i+inbetween)),cy*(sin_scale(i+inbetween)),type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="",xlab="",ylab="",col=c(200,200,200,200))
par(new=T)
}

plot(cx,cy,type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="Sine Scale",xlab="Cos(θ)",ylab="Sin(θ)i", col=rgb(0,1,0))

abline(h=0,col=rgb(1,0,0))
abline(v=0,col=rgb(1,0,0))

# 凡例を書き添える 。
legend("topright", legend=c("Circle(Radius=1,Augment=θ)","x=y=0"), lty =c(1,1),col=c(rgb(0,1,0),rgb(1,0,0)))
}

sin_scale01(0)

アニメーション化(0→1、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/90,2/90,0,1/90,2/90,0,1/90,2/90,0,1/90,2/90)
saveGIF({
for (i in Time_Code){
　sin_scale01(i)
}
}, interval = 0.1, movie.name = "sin_scale01_back.gif")

アニメーション化(1→0、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,2/90,1/90,0,2/90,1/90,0,2/90,1/90,0,2/90,1/90)
saveGIF({
for (i in Time_Code){
　sin_scale01(i)
}
}, interval = 0.1, movie.name = "sin_scale01_forward.gif")

アニメーション化(どちらにも見える、中割り一枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/60,0,1/60,0,1/60,0,1/60,0,1/60,0,1/60)
saveGIF({
for (i in Time_Code){
　sin_scale01(i)
}
}, interval = 0.1, movie.name = "sin_scale01_even.gif")

射影対象をexp(-Sin(θ))と置くとまた別の景色が浮かび上がります。
Rplot19.png

sin_scale <- function(x) Im((0+1i)^(2*x))

sin_scale02<-function(inbetween){

c0<-seq(0,pi*2,length=61)
cx<-cos(c0)
cy<-sin(c0)

c0<-seq(0,1/2,length=31)

for (i in c0){
plot(cx*(exp(-sin_scale(i+inbetween))),cy*(exp(-sin_scale(i+inbetween))),type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="",xlab="",ylab="",col=c(200,200,200,200))
par(new=T)
}

plot(cx,cy,type="l",xlim=c(-1,1),asp=1,ylim=c(-1,1),main="exp(Sine Scale)",xlab="Cos(θ)",ylab="Sin(θ)i", col=rgb(0,1,0))

abline(h=0,col=rgb(1,0,0))
abline(v=0,col=rgb(1,0,0))

# 凡例を書き添える 。
legend("topright", legend=c("Circle(Radius=1,Augment=θ)","x=y=0"), lty =c(1,1),col=c(rgb(0,1,0),rgb(1,0,0)))
}

sin_scale02(0)

アニメーション化(1→exp(-1)、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/180,2/180,0,1/180,2/180,0,1/180,2/180,0,1/180,2/180)
saveGIF({
for (i in Time_Code){
　sin_scale02(i)
}
}, interval = 0.1, movie.name = "sin_scale02_back.gif")

アニメーション化(exp(-1)→1、中割り二枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,2/180,1/180,0,2/180,1/180,0,2/180,1/180,0,2/180,1/180)
saveGIF({
for (i in Time_Code){
　sin_scale02(i)
}
}, interval = 0.1, movie.name = "sin_scale02_forward.gif")

アニメーション化(どちらにも見える、中割り一枚)
image.gif

#アニメーション
library("animation")
Time_Code=c(0,1/120,0,1/120,0,1/120,0,1/120,0,1/120,0,1/120)
saveGIF({
for (i in Time_Code){
　sin_scale02(i)
}
}, interval = 0.1, movie.name = "sin_scale02_even.gif")

考えてみればこれって世界卵(World Egg)の下半分そのものですね。という事は？
【無限遠点を巡る数理】「名状し難きもの(The Unnamable)問題」を世に解き放つ。
image.gif

#相補誤差関数(ELFC=Complementary Error Function)
相補誤差関数(ERFC=Complementary Error Function) erfc(x)=1-erf(x)==2/sqrt(π)×integrate(exp(-t^2), x, ∞) 
image.png

erfc <- function(x) 2 * pnorm(x * sqrt(2), lower=FALSE)
plot(erfc,xlim=c(-3,3),ylim=c(0,2),main="Complementary Error Function(ELFC)",xlab="x",ylab="ELFC(x)",col=rgb(1,0,0))
abline(h=1)
abline(v=0)

x軸とy軸を逆転させた逆関数はこんな感じになります。
image.png

erf <- function(x) 2 * pnorm(x * sqrt(2)) - 1
cy<-seq(-3,3,length=61)
cx<-erfc(cy)
plot(cx,cy,type="l", xlim=c(0,2),ylim=c(-3,3),main="Inverse ELFC",xlab="x",ylab="Inverse ELFC(x)",col=rgb(1,0,0))
abline(h=0)
abline(v=1)

これを同様に円形尺に射影するとこうなります。
【初心者向け】「観測者問題」のあっけない解決方法? - Qiita
image.png

erfc <- function(x) 2 * pnorm(x * sqrt(2), lower=FALSE)

erfc_scale00<-function(inbetween){

c0<-seq(0,pi*2,length=61)
cx<-cos(c0)
cy<-sin(c0)

c0<-seq(0,3,length=31)

for (i in c0){
plot(cx*(1-erfc(i+inbetween))+1,cy*(1-erfc(i+inbetween)),type="l",xlim=c(0,2),asp=1,ylim=c(-1,1),main="",xlab="",ylab="",col=c(200,200,200,200))
par(new=T)
}

plot(cx+1,cy,type="l",xlim=c(0,2),asp=1,ylim=c(-1,1),main="Complementary Error Function(ElfC) Scale",xlab="Cos(θ)",ylab="Sin(θ)i", col=rgb(0,1,0))

abline(h=0,col=rgb(1,0,0))
abline(v=1,col=rgb(1,0,0))

# 凡例を書き添える 。
legend("topright", legend=c("Circle(Radius=1,Augment=θ)","x=1, y=0"), lty =c(1,1),col=c(rgb(0,1,0),rgb(1,0,0)))
}

erfc_scale00(0)

アニメーション化(0→∞、中割り二枚)
image.gif


#アニメーション
library("animation")
Time_Code=c(0,1/30,2/30,0,1/30,2/30,0,1/30,2/30,0,1/30,2/30)
saveGIF({
for (i in Time_Code){
　erfc_scale00(i)
}
}, interval = 0.1, movie.name = "erfc_scale_back.gif")

アニメーション化(∞→0、中割り二枚)
image.gif


#アニメーション
library("animation")
Time_Code=c(0,2/30,1/30,0,2/30,1/30,0,2/30,1/30,0,2/30,1/30)
saveGIF({
for (i in Time_Code){
　erfc_scale00(i)
}
}, interval = 0.1, movie.name = "erfc_scale_forward.gif")

アニメーション化(どちらにも見える、中割り一枚)
image.gif


#アニメーション
library("animation")
Time_Code=c(0,1/20,0,1/20,0,1/20,0,1/20,0,1/20,0,1/20)
saveGIF({
for (i in Time_Code){
　erfc_scale00(i)
}
}, interval = 0.1, movie.name = "erfc_scale_even.gif")

#複素共役(Complex Conjugate)の正面図？
ところで複素共役（複素共軛, Complex Conjugate）の振る舞いは、上のアニメーションの「側面図」とも見て取れるのです。
【初心者向け】複素共役のアニメーション表示について。 - Qiita
image.gif
x軸に沿っての動き。
アニメーション化(x=+1→-1、中割り二枚)
image.gif
アニメーション化(-1→+1、中割り二枚)
image.gif
アニメーション化(どちらにも見える、中割り一枚)
image.gif
y軸に沿っての動き。
アニメーション化(y=+1→-1、中割り二枚)
image.gif
アニメーション化(y=-1→+1、中割り二枚)
image.png
アニメーション化(どちらにも見える、中割り一枚)
image.gif
何かを連想させますね。そう正面図が円描写関数Cos(θ)+Sin(θ)i、側面図がそれぞれCos(θ)とSin(θ)となる単位円筒(Unit Cylinder)概念…
【初心者向け】物理学における「単位円筒」の概念について。
XY軸(円弧)
image.gif
XZ軸(Cos波)
image.gif
YZ軸(Sin波)
image.gif

#そして「統計学」の世界へ…
ちなみに統計学ではしばしば「95%棄却域」や「99%棄却域」なる概念を用いますが、それに対応するのはおそらく以下です。
23-5. 検定統計量と棄却域・採択域 | 統計学の時間 | 統計WEB

#「95%棄却」ではelf(x)のxの値が概ねsqrt(2)に該当。

sqrt(2)
[1] 1.414214

#elf(x)=0.95におけるxの下限は1.3859~1.3860

> erf(1.3859)
[1] 0.9499994
> erf(1.3860)
[1] 0.9500159

#elf(x)=0.95におけるxの上限は1.4522~1.4523

> erf(1.4522)
[1] 0.9599973

 > erf(1.4523)
[1] 0.960011

#elf(x)=0.99におけるxの下限は1.8213~1.8214

> erf(1.8213)
[1] 0.9899965

> erf(1.8214)
[1] 0.9900006

確かに恐ろしく収束が早く「誤差として切って捨てる」のが惜しくありません。一方、後世の統計学はそれを「決して偶然起こった事ではない」基準として採用した訳です。

#Pythonの場合

import math as m
import cmath as c
import numpy as np
import scipy.special as sp
import matplotlib.pyplot as plt
import matplotlib.animation as animation

#単位円データ作成
c0=np.linspace(0,2*m.pi,61,endpoint = True)
s0=[]
for nm in range(len(c0)):
    s0.append(complex(m.cos(c0[nm]),m.sin(c0[nm])))
s1=np.array(s0)

#複素共役データ作成

#描画準備
plt.style.use('default')
fig = plt.figure(111)

#関数定義
def sin_scale(x):
    return ((0+1j)**(2*x)).imag
def erf_scale00(n):
    plt.cla()
    #共益線描画
    z0=np.linspace(0,3,31)
    for nm in z0:
        plt.plot(s1.real*sp.erf(nm+Time_code[n]),s1.imag*sp.erf(nm+Time_code[n]),color="black",lw=0.5);    
    #円周描画
    plt.plot(s1.real,s1.imag,color="green", label="Unit Cylinder")

    plt.ylim([-1.1,1.1])
    plt.xlim([-1.1,1.1])
    plt.title("Error Function(Elf)")
    plt.xlabel("Real")
    plt.ylabel("Imaginal")
    ax = fig.add_subplot()
    ax.set_aspect('equal', adjustable='box')
    ax.legend(loc='upper right')
    #補助線描画
    plt.axvline(0, 0, 1,color="red")
    plt.axhline(0, 0, 1,color="red")
    
#アニメーション化(0→∞、中割り二枚)
Time_code=[0,1/30,2/30,0,1/30,2/30,0,1/30,2/30,0,1/30,2/30]
#アニメーション化(∞→0、中割り二枚)
#Time_code=[0,2/30,1/30,0,2/30,1/30,0,2/30,1/30,0,2/30,1/30]
#アニメーション化(どちらにも見える、中割り一枚)
#Time_code=[0,1/20,0,1/20,0,1/20,0,1/20,0,1/20,0,1/20]
#erf_scale00(1)
#plt.show()

ani = animation.FuncAnimation(fig, erf_scale00, interval=50,frames=len(Time_code))
ani.save("erf008.gif", writer="pillow")

アニメーション化(0→∞、中割り二枚)
image.gif
アニメーション化(∞→0、中割り二枚)
image.gif
アニメーション化(どちらにも見える、中割り一枚)
image.gif

そんな感じで以下続報…

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?