統計で使える身近な関数をPythonで作ってみた

##目次
- この記事を書こうと思ったきっかけ
- 実行環境
- 平均
- レンジ
- 中央値1
- 中央値2(昇順、降順の時のみ)
- 標準偏差
- 分散
- 累積相対度数
- 変動係数
- 2変量解析
    - 共分散
    - 相関係数
- おまけ編

##この記事を書こうと思ったきっかけ
図書館で統計学の本を読んでいたところ、「Pythonを使って(標準偏差などを)求められるのでは」と思ったこと。また、Pythonの実力をあげたかったので、気分転換も兼ねてこの記事を書きました。


##実行環境
macOS Big Sur
11.1
paiza.IO(Python3)で実行しています。


##平均
####プログラム
```Python
def mean(l):
    return sum(l)/len(l)
```
####使用例
```Python
print(mean([3,4,5,7,6]))
#5.0
print(mean([3,4,11,7,6]))
#6.2
```
####注意点
- 戻り値はfloat型
- 昇順でないリストにも使える点
- 数がとてつもなく大きくなると下のように「e」で省略されてしまう点

```Python
print(mean([3,4,1111111111111111111111111111111,7,6]))
#2.222222222222222e+29
```

##レンジ
####プログラム
```Python
def Range(l):
    return max(l)-min(l)
```
昇順でなくても良い。

####注意点
- 組込み変数のrangeと重ならないように関数名の1文字目を大文字にしています。

##中央値1
####プログラム
```Python
def chuuouti(l):
    l2 = sorted(l)
    if len(l2)%2:
        return l2[int(len(l)/2)]
    else:
        return (l2[int(len(l)/2)-1]+l2[int(len(l)/2)])/2
```
####使用例
```Python
print(chuuouti([1,4,2,3]))
```

####注意点
- 要素数が偶数の時は真ん中の2つの要素の平均をとる。
- 昇順でなくても良い(しかし、sortedの処理に少し時間がかかるので、昇順の時は「中央値2」を使うと良い。

####例題
[第四回 アルゴリズム実技検定　A - 中央値](https://atcoder.jp/contests/past202010-open/tasks/past202010_a)

##中央値2(昇順、降順の時のみ)
####プログラム
```Python
def chuuouti(l):
    if len(l)%2:
        return l[int(len(l)/2)]
    else:
        return (l[int(len(l)/2)-1]+l[int(len(l)/2)])/2
```
####注意点
- 昇順、降順の時以外は使えない
- 長さが偶数である時は戻り値がfloat型だが、長さが奇数の時は戻り値の型が中央値と一致する点。(下)

```Python
#長さが偶数の時
print(chuuouti([3,2,2,1]))
#2.0
print(chuuouti([4,3,2,1]))
#2.5

#長さが奇数の時
print(chuuouti([1,3,7]))
#3
print(chuuouti([1,3.2,7]))
#3.2
```


##標準偏差
####プログラム
```Python
from math import sqrt
def deviation_value(l):
    mean = sum(l)/len(l)#平均
    deviation_square = [(l[i]-mean)**2 for i in range(len(l))]#偏差平方
    deviation_square_sum = sum(deviation_square)#偏差平方和
    return sqrt(deviation_square_sum/len(l))#偏差値
```
####使用例
```Python
print(deviation_value([10,15,5,-22,-18]))
#15.086417732516887
```
上の使用例では平均が{10+15+5+(-22)+(-18)}/5より-2.
よって、偏差平方は144,289,49,400,256.この和は1138.
これを要素数(5個)で割り、平方根をとれば、標準偏差になる。

##累積相対度数
####プログラム
```Python
def cumulative_relative_frequency(x):
    wa = sum(x)#全体の和
    l = [0]
    y = [x[i]/wa for i in range(len(x))]
    #ここから累積和
    for i in range(len(y)):
        l.append(l[-1]+y[i])
    del l[0]
    return l
```

####使用例
```Python
print(cumulative_relative_frequency([1,3,7,6,9]))
#[0.038461538461538464, 0.15384615384615385, 0.4230769230769231, 0.6538461538461539, 1.0]
```

##変動係数
####プログラム
```Python
from math import sqrt
def coefficient_of_variation(l):
    mean = sum(l)/len(l)#平均
    deviation_square = [(l[i]-mean)**2 for i in range(len(l))]#偏差平方
    deviation_square_sum = sum(deviation_square)#偏差平方和
    return sqrt(deviation_square_sum/len(l))/mean
```
####使用例
```Python
print(coefficient_of_variation([4,5,3,5]))
#https://bellcurve.jp/statistics/course/19515.htmlの例。
```
本当は0.2であるが、0.195と出力しているため、少し誤差があるとみられる。この誤差は、平方根(sqrt)の時に発生したものとみられる。

#2変量解析
##共分散
####プログラム

```Python
def bivariate_analysis(x,y):
    mean_x = sum(x)/len(x)
    mean_y = sum(y)/len(y)
    deviation_x = [x[i]-mean_x for i in range(len(x))]
    deviation_y = [y[i]-mean_y for i in range(len(y))]
    product_of_deviations = [deviation_x[i]*deviation_y[i] for i in range(len(x))]
    return sum(product_of_deviations)/len(x)
```
####使用例
```Python
print(bivariate_analysis([4,7],[3,6]))
```
この関数を使うことで、2つのリストの共分散が求められる。

####メリット
- 正確な値を求められる。

##相関係数
####プログラム
```Python
def correlation_coefficient(x,y):
    return bivariate_analysis(x,y)/deviation_value(x)/deviation_value(y)
```
####標準偏差と共分散を入れたプログラム
```Python
from math import sqrt

def deviation_value(l):#標準偏差
    mean = sum(l)/len(l)#平均
    deviation_square = [(l[i]-mean)**2 for i in range(len(l))]#偏差平方
    deviation_square_sum = sum(deviation_square)#偏差平方和
    return sqrt(deviation_square_sum/len(l))#標準偏差
    
    
def bivariate_analysis(x,y):#共分散
    mean_x = sum(x)/len(x)
    mean_y = sum(y)/len(y)
    deviation_x = [x[i]-mean_x for i in range(len(x))]
    deviation_y = [y[i]-mean_y for i in range(len(y))]
    product_of_deviations = [deviation_x[i]*deviation_y[i] for i in range(len(x))]
    return sum(product_of_deviations)/len(x)
    
def correlation_coefficient(x,y):
    return bivariate_analysis(x,y)/deviation_value(x)/deviation_value(y)
```
####注意点
- 標準偏差と共分散のプログラムを使っています。

##おまけ編
n人のクラスがあります。誕生日が同じの生徒が存在する確率を求めよ。

##<おまけ編の実装編>
[同じ誕生日の二人組がいる確率について](https://mathtrain.jp/birthday)をつかって実装していきたいと思います。[同じ誕生日の二人組がいる確率について](https://mathtrain.jp/birthday)より、求める確率は

```math
1−_{365}P_n\div365^n
```
になる。
#####プログラム
```Python
from math import factorial
def permutations_count(n, r):#下の※のリンクから抜粋
    return factorial(n) // factorial(n - r)
def birthday(n):#n人のクラスで同じ誕生日の人がいる確率
    return 1-permutations_count(365,n)/365**n
```

※[Pythonで階乗、順列・組み合わせを計算、生成](https://note.nkmk.me/python-math-factorial-permutations-combinations/)からP(n,r)のコードを抜粋した。

このプログラムを使って初めて75%を超えるのはnがいくつ以上の時であるか、求めてみましょう。また、その時の確率を求めてみましょう。

####プログラム
```Python
people = 0
while True:
    if birthday(people) >= 0.75:
        print(people)
        print(birthday(people))
        break
    people += 1
```
peopleはクラスの人数です。birthday(100)が75%を超えると言うことがわかっていれば、for文でも良いかもしれません。23人ですでに50.7%のようです。




##参考文献
- 翻訳
    - [goo辞書](https://dictionary.goo.ne.jp)
    - [Google翻訳](https://translate.google.com/?hl=ja&sl=ja&tl=en&op=translate)
- 書き方
    - [色々な確率分布の起こる現象を理解した上でpythonで実装
](https://qiita.com/phyblas/items/226132b3c1ec5467b267)
- プログラム
    - [
The Programing Guide]
(https://strawberryprogrami.wixsite.com/programing)
- 用語の意味
    - [すぐに使える統計学　菅 民郎,土方 裕子著　初版第1刷 出版社:SoftBank Creative](https://www.amazon.co.jp/すぐに使える統計学-菅-民郎/dp/4797349077)
    - [共分散の意味と簡単な求め方](https://mathtrain.jp/covariance)
    - [5-3. 変動係数を求めてみよう](https://bellcurve.jp/statistics/course/19515.html)
    - [数学記号（等号、不等号、演算子、集合）](https://medemanabu.net/latex/operators/)
    - 図解雑学　統計　今野紀雄著　ナツメ社