乙Py先生のプログラミング教室

乙Py先生のプログラミング教室
初学者のためのプログラミング学習サイト

今回は日本で一番人気の漫画である
Hunter x Hunterに登場するアイテム
「リスキーダイス」をMatplotlibで
制作してみました

解説動画はこちら



リスキーダイスについての説明は
改めてする必要も無いでしょうけども
一応やっておきますね

HUNTER × HUNTER に出てくる
アイテムの一つで20面体のサイコロ

19面が大吉で、1面が大凶
大吉が出るととてもいいことが起こる
ただし大凶が出るとそれまでに出た
大吉分がチャラになるほどの不幸が起きる

という素敵なアイテムです

これを作っていくのを考えましょう

正20面体は正三角形が
20個組み合わさった多面体です

これを展開したものを
Matplotlibで描写していきます

正三角形を描くには
3点が必要です

三角形

頂点の座標
正三角形の高さは
公式で求められます

それぞれの座標を求めて3角形を描き
展開図になるように3角形を
配置すれば良いですね

ただそれだけだと組み立てられないので
余白、のり代部分を作ってあげます

Matplotlibで描画する際には
三角形を描くメソッド
matplotlib.patches.Polygon
というのが有ります

これに3点を渡してあげれば
三角形が描画出来ます

あとは「大吉」などの文字を
ax.text で描画してあげれば
20面体のサイコロとなりますが
Matplotlibは初期状態だと
日本語フォントに対応してないので
別途インストールが必要です


リスキーダイスの展開図を描画するコードはこちら
import matplotlib.pyplot as plt
import matplotlib.patches as pat

fig = plt.figure(figsize=(16, 16))
ax = plt.axes()

# 正20面体
l=1
h = (l**2 - (l/2)**2) ** 0.5
for i in range(5):
    x=l*i
    # 1段目
    p = pat.Polygon(xy = [(x,0), (x+1,0), (x+l/2 , h)],fc = "gray",ec = "black")
    ax.add_patch(p)
    p = pat.Polygon(xy = [(x+l/2,h), (x+l/2+l,h), (x+l , 0)],fc = "white",ec = "black")
    ax.add_patch(p)
    # 2段目
    p = pat.Polygon(xy = [(x+l/2,h), (x+l/2+l,h), (x+l , h*2)],fc = "gray",ec = "black")
    ax.add_patch(p)
    p = pat.Polygon(xy = [(x,h*2), (x+l,h*2), (x+l/2 , h)],fc = "gray",ec = "black")
    ax.add_patch(p)
    # 3段目
    p = pat.Polygon(xy = [(x,h*2), (x+1,h*2), (x+l/2 , h*3)],fc = "white",ec = "black")
    ax.add_patch(p)
    p = pat.Polygon(xy = [(x+l/2,h*3), (x+l/2+l,h*3), (x+l , h*2)],fc = "gray",ec = "black")
    ax.add_patch(p)

# のり代
p = pat.Polygon(xy = [(l*5,h*2), (l*5+l,h*2), (l*5+l/2 , h)],fc = "white",ec = "black")
ax.add_patch(p)

# 切り取りライン
for i in range(4):
    x=l*i + l
    plt.plot([x, x+l/2], [0, h],color="red")
    plt.plot([x, x+l/2], [h*2, h*3],color="red")

# テキスト
for i in range(5):
    x=l*i
    if i!=0:
        ax.text(x + l/2 , h/3 , '大吉',color="white",size=30,va='center', ha='center')
    else:
        ax.text(x + l/2 , h/3 , '大凶',color="red",size=30,va='center', ha='center')
    ax.text(x + l , h/3+h , '大吉',color="white",size=30,va='center', ha='center')
    
    # 逆さ文字
    ax.text(x + l , h*3-h/3 , '大吉',color="white",size=30,va='center', ha='center',rotation = 180)
    ax.text(x + l/2 , h*2-h/3 , '大吉',color="white",size=30,va='center', ha='center',rotation = 180)

plt.axis('scaled')
plt.axis('off')
plt.show()

実行すると・・・
こんな感じのが出来上がります

リスキーダイス

リスキーダイス


これをカラー印刷すると・・・
リスキーダイス展開図

いらない部分を切り取って
赤い線は切り取り線で
白い部分がのり代です

いい感じに糊付けすると・・・

リスキーダイス完成形

いい感じにサイコロになりました

W杯
日本 vs クロアチア

何度か振った結果・・・

大吉しか出ませんでしたーーー


なので、クロアチアにも
勝ってくれる事でしょうwww

今回はリスキーダイスを
Matplotlibで制作してみました

作りたい方は画像をダウンロードして
印刷して工作してみて下さい

それでは
 

今回はプロジェクトオイラーの
数学の問題を解いてみました

解説動画はこちら


 
さて
プロジェクト・オイラー

は、英語で書かれた
数学的/コンピューター プログラミングの問題集です
リンクから色々な問題をみる事ができます

レオンハルト・オイラーは
18世紀の数学者・天文学者のことですね

今回はこれを解いていきたいと思います
数学の問題をプログラムを
駆使して解いてみましょう

問題を考えたい人は
動画を止めて考えてみて下さい

解答そのものは動画をご覧ください






問題1

1000未満の3か5の倍数の合計値を求めよ






・・・・・

・・・・

・・・

・・










解答1

比較的簡単な問題です

1000未満の数値の生成
3と5の倍数の判定
合計値を求める

これらが達成出来ていれば
解けると思います

sum([i for i in range(1,1000) if i%3==0 or i%5==0])





問題2

400万未満のフィボナッチ数で

偶数の物の合計値を求めよ

フィボナッチ数とは

1, 2, 3, 5, 8, 13, 21, 34, 55, 89 のように

どの数字も前2つの数字を足した数字のこと




・・・・・

・・・・

・・・

・・



解答2

フィボナッチ数を求める関数を作って
400万までのリストを作り
偶数のものを抜き出せば解けます

フィボナッチ数のリストを返す関数
# その数までのフィボナッチ数を求める
def fib(num):
    a , b , fib_l = 0 , 1 , []
    while b <= num:
        fib_l.append(b)
        a , b = b , a+b
    return fib_l

sum(([f for f in fib(4000000) if f%2==0]))

ジェネレーター式でも同様に出来ます
# ジェネレーター形式
def fib_gen(num):
    a , b = 1 , 1
    while b <= num:
        yield b
        a , b = b , a+b
        
sum(([f for f in fib_gen(4000000) if f%2==0]))




問題3

600851475143 の最大の素因数を求めよ

素因数分解をして、その最大値を求める


・・・・・

・・・・

・・・

・・






解答3

素因数分解をする手法はたくさん有りますが
中でも簡単な試し割り法を使って
素因数分解してみましょう

素因数分解しようとする整数nを
小さい順に割ってみて
割り切れるかどうかを調べる手法です


num = 600851475143

# 3から奇数で割っていく
f , a = 3 , []
while f * f <= num:
    
    # 割り切れたら素因数に採用
    if num % f == 0:
        a.append(f)
        
        # 次の数は、素因数で割った数で探す
        num //= f
    
    # 割り切れなければ次の素数へ
    else:
        f += 2
        
if num != 1:
    a.append(num)

# 素因数の組み合わせ
print(a)

# 素因数の最大値 
print(sorted(a,reverse=True)[0])


なお、ライブラリを用いれば
簡単に素因数分解できてしまいます
# ライブラリで素因数分解
import sympy

num = 600851475143
sympy.factorint(num)





問題4

2 つの 3 桁の数の積から作られる

最大の回文数値を求めよ

101*101 = 10201 のような

回文になる数値の最大値を求める


・・・・・

・・・・

・・・

・・






解答4

数値を文字列に直してひっくり返せば
回文の出来上がりです

answer = {}

for a in range(100,1000):
    for b in range(100,1000):
        s_num = str(a * b)
        if s_num == s_num[::-1]:
            answer[a * b] = [a,b]

for k,v in sorted(answer.items(),reverse=True):
    print(k,v)
    break








問題5

1 から 20 までのすべての数で割り切れる

最小の正の数を求めよ

最小公倍数を求める問題


・・・・・

・・・・

・・・

・・







解答5

最小公倍数は
次の様な計算式で求められます

最小公倍数 = 数値の積 ÷ 最大公約数


最大公約数を求めるには
ユークリッドの互除法を用います

ユークリッドの互除法は
2つの自然数の最大公約数を求める
手法の一つですが
比較的実装が楽です

# ユークリッドの互除法
# 最大公約数を求める
def gcd(a, b):
    if a == 0:
        return b
    m = b % a
    return gcd(m, a)

# 最小公倍数 = 数値の積 / 最大公約数
from functools import reduce

# 最小公倍数を計算する(2つ)
def lcm_base(x, y):
    return (x * y) // gcd(x, y)

# 最小公倍数を計算する(複数)
def lcm(*nums):
    return reduce(lcm_base , nums, 1)


s = [i for i in range(2,21)]
# 最小公倍数をリストで求める際は * を付ける
lcm(*s)



python3.8以降なら
ライブラリでも求められます
# 最小公倍数を求める
import math
s = [i for i in range(2,21)]

# 最小公倍数をリストで求める際は * を付ける
math.lcm(*s)






問題6

100 個の自然数の

和の 2 乗と 2 乗の和の差を求めよ

例1,2,3

(1 + 2 + 3)^2 − (1^2 + 2^2 + 3^2) 

= 36 −14 

= 22



・・・・・

・・・・

・・・

・・






解答6


和の2乗は比較的楽に計算できますが
2乗の和は少しコツがいります

map関数は関数を引数に指定できるので
リストの要素に
関数を適応した結果を求めるのに
すごく役立ちます
nums = [i for i in range(1,101)]

# 和の2乗
sum_square = sum(nums) ** 2
print(sum_square)

# 2乗の和
sum_of_squares = sum(map(lambda n : n**2 , nums))
print(sum_of_squares)

# 2乗の和と和の2乗の差
print(sum_square - sum_of_squares)





問題7

10001 番目の素数を求めよ





・・・・・

・・・・

・・・

・・





解答7

素数を求めるのは流石に面倒くさいので
ライブラリの力を借りましょう

isprimeメソッドは
実行速度がめちゃくちゃ速いので
大抵の人が実装したものよりも
速いだろうと思います
# 素数を求めるライブラリ
from sympy import isprime

primes = [2,3]
for i in range(5,200001,2):
    if isprime(i):
        primes.append(i)
    if len(primes)==10001:
        answer = i
        print(answer)
        break


ライブラリでその数値までの
素数の個数を数える事も出来ます
# ライブラリで素数の個数を数える
from sympy import primepi

print(primepi(answer))







問題8

次の1000 桁の数の中で

最大の積を持つ隣接する13桁は何か

数値は下記をお使いください

73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450






・・・・・

・・・・

・・・

・・





解答8


まずは文字列をデータ化しますが
きちんとゴミを取り除いておきましょう

文字列を数値に直して
計算すればOKです

総積を求めるのは
for文を使っても良いですが
reduce関数などを使って
まとめる事ができます


data = '''
73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450
'''.replace('\n','')

import operator
from functools import reduce

answer = {}
for i in range(0 , len(data)-12):
    nums = data[i : i+13]
    if '0' in nums:
        continue
    
    # 総積を求める
    a = reduce(operator.mul, [int(n) for n in nums],1)
    answer[a] = nums
    
for k,v in sorted(answer.items(),reverse=True):
    print('数値 :  {0} , 文字列 : {1}'.format(k,v))
    break







問題9

a < b < c で a + b + c = 1000である時

ピタゴラスのトリプレットとなる数値の

 a * b * c の値を求める

ピタゴラスのトリプレットとは

a^2 + b^2 = c^2






・・・・・

・・・・

・・・

・・




解答9

a<b<c かつ a+b+c=1000 なので
abcの数値の最大値はおのずと決まってきます

あとはピタゴラスの定理のような
関係性になる数値を探し当てれば
それが答えです


# aの最大値は1000の3分の1より小さい
for a in range(1, 332):
    
    # bの最小値は a+1、最大値は 1000から a を引いた半分
    for b in range(a + 1 , (1000 - a)//2):
        c = 1000 - (a + b)
        if a**2 + b**2 == c**2:
            print(a * b * c , a , b , c)
            break
    else:
        continue
    break





さて、ここまでが今回の内容です
なかなか面白い問題ばかりですねえ


まだまだプロジェクトオイラーには
たくさんの問題が残っているので
これからも問題を解いていきたいと思います

それでは

今回はイロレーティングを使って
日本がW杯で優勝する確率を
シミュレーションしてみました


解説動画はこちら



さて、今回はワールドカップ2022目前ということで
イロレーティングを用いて結果がどうなるのか
シミュレーションしてみることとしました

イロレーティングというのは
対戦ゲームやスポーツに
良く用いられているもので
相対評価で実力を表すために
使われる指標の一つです

レーティングを使うと
対戦国同志の勝率を出すことができます

レーティングから勝率を求める計算式は
次のようになっており

download


対戦国AとBでの
勝つ確率を求める事ができます

レーティングが高い方が
勝つ確率が上がります

おおよそのレーティング差での
勝率は次のようになっています

レート差勝率
1051.44%
5057.15%
10064.01%
20075.97%
30084.9%
40090.91%
50094.68%



今回は2022出場国のレーティングを
掲載しているサイトがあったので
そちらのデータを使うこととします。


レーティングの元はこちら


・データ化

このリンクを元にデータ化します
data="""
グループA	カタール	50	1439.89
グループA	エクアドル	44	1464.39
グループA	セネガル	18	1584.38
グループA	オランダ	8	1694.51
グループB	イングランド	5	1728.47
グループB	イラン	20	1564.61
グループB	アメリカ	16	1627.48
グループB	ウェールズ	19	1569.82
グループC	アルゼンチン	3	1773.88
グループC	サウジアラビア	51	1437.78
グループC	メキシコ	13	1644.89
グループC	ポーランド	26	1548.59
グループD	フランス	4	1759.78
グループD	オーストラリア	38	1488.72
グループD	デンマーク	10	1666.57
グループD	チュニジア	30	1507.54
グループE	スペイン	7	1715.22
グループE	コスタリカ	31	1503.59
グループE	ドイツ	11	1650.21
グループE	日本	24	1559.54
グループF	ベルギー	2	1816.71
グループF	カナダ	41	1475
グループF	モロッコ	22	1563.5
グループF	クロアチア	12	1645.64
グループG	ブラジル	1	1841.3
グループG	セルビア	21	1563.62
グループG	スイス	15	1635.92
グループG	カメルーン	43	1471.44
グループH	ポルトガル	9	1676.56
グループH	ガーナ	61	1393
グループH	ウルグアイ	14	1638.71
グループH	韓国	28	1530.3
""".split('\n')

# テキストからデータ化
data2 = [ d.split('\t') for d in data[1:-1]]
data2 = [d[0:3] + [float(d[3])] for d in data2]

# 国別レーティングのデータ化
data3 = {d[1]:d[3] for d in data2}

レーティングから勝率を計算する
関数はこのようになります
# レーティングから勝率計算の関数
def win_rate(a,b):
    return 1/(10**((b-a)/400) +1)


・予選の予測


さてここからは予選結果を
予測していくこととしましょう

グループリーグは総当たりで
2位までが決勝に残る事ができます

本来は勝ち点による
引き分けや得点差を加味したものなのですが
今回はシンプルにグループリーグの
勝ち残りを選ぶという仕様で
勝ち負けだけで選びます

前提条件:勝ち負けのみ、引き分けなし


さっそくグループを作って
総当たりの組み合わせを見てみましょう

グループ4チームから
総当たりの組み合わせを作るには
itertoolsのcombinationsで出来ます


import itertools
import numpy as np

# グループ分け
groups = [data2[0+i*4 : 4+i*4] for i in range(8)]

# 1グループ選ぶ
group = groups[0]
group_name = group[0][0]

# 総当たり戦の組み合わせを作る
brute_forces = list(itertools.combinations([g[1] for g in group],2))
print(brute_forces)
[('カタール', 'エクアドル'), 
('カタール', 'セネガル'), 
('カタール', 'オランダ'), 
('エクアドル', 'セネガル'), 
('エクアドル', 'オランダ'), 
('セネガル', 'オランダ')]

あとはこのチーム同士のレーティングを取ってきて
勝率にあてはめて、勝敗を決めれば
予選の勝利数が計算できます

その上位2チームが残るという訳です
# リーグ戦の結果
league_result = {g[1]:0 for g in group}
for a,b in brute_forces:
    # 国のレーティングをデータから取る
    rating_a = data3[a]
    rating_b = data3[b]
    
    # レーティングから勝率計算
    rate = win_rate(rating_a,rating_b)

    # 勝率で勝者を選ぶ
    choice = np.random.choice([a,b], p=[rate , 1-rate])
    league_result[choice] +=1
    print(a,b,rate,choice)
    
print(league_result)

wins = sorted(league_result.items(),reverse=True,key=lambda x:x[1])[0:2]
print(group_name,wins)
{'カタール': 2, 'エクアドル': 1, 'セネガル': 2, 'オランダ': 1}
グループA [('カタール', 2), ('セネガル', 2)]

こんな感じで、予選を勝つチームが
計算できました

ここまでを関数化してしまいます
def choice_wins(a,b):
    rating_a = data3[a]
    rating_b = data3[b]
    rate = win_rate(rating_a,rating_b)
    return np.random.choice([a,b], p=[rate , 1-rate])

def choice_group_wins(group):
    brute_forces = list(itertools.combinations([g[1] for g in group],2))
    league_result = {g[1]:0 for g in group}
    for a,b in brute_forces:
        rating_a = data3[a]
        rating_b = data3[b]
        choice = choice_wins(a,b)
        league_result[choice] +=1
    return sorted(league_result.items(),reverse=True,key=lambda x:x[1])[0:2]

# 予選リーグの結果
league_results = []
for group in groups:
    wins = choice_group_wins(group)
    print(wins)
    league_results.append(wins)
[('カタール', 2), ('セネガル', 2)]
[('イングランド', 2), ('アメリカ', 2)]
[('アルゼンチン', 3), ('メキシコ', 2)]
[('デンマーク', 3), ('フランス', 2)]
[('スペイン', 3), ('日本', 2)]
[('ベルギー', 2), ('モロッコ', 2)]
[('ブラジル', 2), ('スイス', 2)]
[('韓国', 3), ('ポルトガル', 1)]


勝率から計算しているので
実行毎にランダムで結果が変わります
これで予選を勝ったチームが出来ました



・予選リーグの突破確率

予選リーグの計算を10000回試行して
日本が予選を突破できる確率を
求めてみましょう

グループEのレーティングは
次の様になっています


スペイン 1715.22
ドイツ 1650.21
日本    1559.54
コスタリカ 1503.59


日本の対戦国との勝率は
次のようになります
# 日本の対戦国との勝率
rate1 = win_rate(data3['日本'],data3['スペイン'])
rate2 = win_rate(data3['日本'],data3['ドイツ'])
rate3 = win_rate(data3['日本'],data3['コスタリカ'])

print('{0}戦の勝率 : {1:.03f}%'.format('スペイン_' , rate1*100))
print('{0}戦の勝率 : {1:.03f}%'.format('_ドイツ_' , rate2*100))
print('{0}戦の勝率 : {1:.03f}%'.format('コスタリカ' , rate3*100))
スペイン_戦の勝率 : 28.984%
_ドイツ_戦の勝率 : 37.240%
コスタリカ戦の勝率 : 57.983%


日本の予選リーグ突破は
次のコードで計算しています
# 予選突破回数
japan_results = {g[1]:0 for g in groups[4]}

for i in range(10000):
    wins = choice_group_wins(groups[4])
    japan_results[wins[0][0]]+=1
    japan_results[wins[1][0]]+=1

for k,v in sorted(japan_results.items(),reverse=True,key=lambda x:x[1]):
    print('レーティング : {0:.02f} 、突破回数 : {1:04} , {2}'.format(data3[k] , v , k))
レーティング : 1715.22 、突破回数 : 8267 , スペイン
レーティング : 1650.21 、突破回数 : 5877 , ドイツ
レーティング : 1503.59 、突破回数 : 3030 , コスタリカ
レーティング : 1559.54 、突破回数 : 2826 , 日本


毎回変わりますが
おおよそ30%ほどの確率で
予選突破できるのではないかと思われます




・決勝の計算

決勝の組み合わせは次のようになっていて
各グループの1位と2位が戦う仕組みです

download-1

これをプログラムに起こしていくと
次のようになります
def final_wins(league_results):
    # 決勝リーグ1回戦の組み合わせ
    final_stage_1_result = []
    for i in range(0,7,2):
        a1 = league_results[i][0][0]
        b2 = league_results[i+1][1][0]
        choice1= choice_wins(a1,b2)
        final_stage_1_result.append(choice1)

        a2 = league_results[i][1][0]
        b1 = league_results[i+1][0][0]
        choice2 = choice_wins(a2,b1)
        final_stage_1_result.append(choice2)

    # 準々決勝の組み合わせ
    quater_finals = []
    for i in [0,4]:
        a1 = final_stage_1_result[i]
        b2 = final_stage_1_result[i+2]
        choice1= choice_wins(a1,b2)
        quater_finals.append(choice1)

        a2 = final_stage_1_result[i+1]
        b1 = final_stage_1_result[i+3]
        choice2 = choice_wins(a2,b1)
        quater_finals.append(choice2)
        
    # 準決勝の組み合わせ
    semi_finals = []
    a1 = quater_finals[0]
    b2 = quater_finals[2]
    choice1= choice_wins(a1,b2)
    semi_finals.append(choice1)

    a2 = quater_finals[1]
    b1 = quater_finals[3]
    choice2 = choice_wins(a2,b1)
    semi_finals.append(choice2)
    
    # 決勝の組み合わせ
    a1 = semi_finals[0]
    b2 = semi_finals[1]
    return choice_wins(a1,b2)

これで優勝国が導き出されます

これを使って優勝確率を
計算していきましょう


・日本が優勝できる確率は?


ここまでを10000回試行して


final_calc = {k:0 for k in data3.keys()}
for i in range(10000):
    # 予選リーグの結果
    league_results = []
    for group in groups:
        wins = choice_group_wins(group)
        league_results.append(wins)
    # 優勝国
    winner = final_wins(league_results)    
    final_calc[winner] +=1

for k,v in sorted(final_calc.items(),reverse=True,key=lambda x:x[1]):
    print('レーティング : {0:.02f} 、優勝回数 : {1:04} , {2}'.format(data3[k] , v , k))


結果はプログラムを動かすか
動画をご覧ください


まあ、普通にレーティングが高い国が
優勝する確率も高い、という
当たり前の結果になります

今回はイロレーティングを用いて
W杯で優勝する確率を
求めてみました

それでは。




 

このページのトップヘ