EM 알고리즘

import numpy as np
import pandas as pd

data = pd.read_csv('./em_data.csv').drop('Unnamed: 0',axis=1).fillna(3)

def get_count(i, j, k, data):
    a = data.loc[((data.y1==i) & (data.y2==j) & (data.y3==k)),'Count']
    return float(a)

pi = np.random.rand(8).reshape(2,2,2)
pi = pi/np.sum(pi)
temp_count = pi
ijk = [(1, 1, 1), (1, 1, 0), (1, 0, 1), (1, 0, 0),
       (0, 1, 1), (0, 1, 0), (0, 0, 1), (0, 0, 0)]

error = 10
cnt = 0
while error>1e-10:
    pi_old = pi
    cnt = cnt+1
    # E-step and M-step
    for i, j, k in ijk:
        temp_count[i, j, k] = get_count(i, j, k, data) +\
            get_count(i, 3, k, data)*pi[i, j, k]/(pi[i, 0, k]+pi[i, 1, k]) +\
            get_count(3, j, k, data)*pi[i, j, k]/(pi[0, j, k]+pi[1, j, k])
    pi = temp_count / np.sum(temp_count)
    # calculate error
    error = np.sum(np.abs(pi_old - pi))
    # print result
    print('-------{}times-------'.format(cnt))
    for i, j, k in ijk:
        print('Pi_({},{},{}) : {:.3}'.format(i,j,k,pi[i,j,k]))
    print()
print('EM iteration {}times.'.format(cnt))

-------1times-------
Pi_(1,1,1) : 0.48
Pi_(1,1,0) : 0.00684
Pi_(1,0,1) : 0.201
Pi_(1,0,0) : 0.00839
Pi_(0,1,1) : 0.0604
Pi_(0,1,0) : 0.0057
Pi_(0,0,1) : 0.218
Pi_(0,0,0) : 0.0194

-------2times-------
Pi_(1,1,1) : 0.381
Pi_(1,1,0) : 0.00772
Pi_(1,0,1) : 0.239
Pi_(1,0,0) : 0.00755
Pi_(0,1,1) : 0.0437
Pi_(0,1,0) : 0.00485
Pi_(0,0,1) : 0.292
Pi_(0,0,0) : 0.0246

-------3times-------
Pi_(1,1,1) : 0.374
Pi_(1,1,0) : 0.00825
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00672
Pi_(0,1,1) : 0.0369
Pi_(0,1,0) : 0.00423
Pi_(0,0,1) : 0.302
Pi_(0,0,0) : 0.0255

-------4times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00865
Pi_(1,0,1) : 0.243
Pi_(1,0,0) : 0.0063
Pi_(0,1,1) : 0.0344
Pi_(0,1,0) : 0.00391
Pi_(0,0,1) : 0.304
Pi_(0,0,0) : 0.0258

-------5times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00889
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00608
Pi_(0,1,1) : 0.0335
Pi_(0,1,0) : 0.00374
Pi_(0,0,1) : 0.304
Pi_(0,0,0) : 0.026

-------6times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00902
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00595
Pi_(0,1,1) : 0.0331
Pi_(0,1,0) : 0.00364
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0261

-------7times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00909
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00589
Pi_(0,1,1) : 0.033
Pi_(0,1,0) : 0.00359
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0261

-------8times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00913
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00585
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00356
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------9times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00915
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00583
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00354
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------10times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00917
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00582
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00353
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------11times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00917
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00581
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00353
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------12times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00581
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00353
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------13times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00581
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00353
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------14times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.00581
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------15times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------16times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------17times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------18times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------19times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------20times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------21times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------22times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------23times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------24times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------25times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------26times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------27times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------28times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------29times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------30times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------31times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

-------32times-------
Pi_(1,1,1) : 0.375
Pi_(1,1,0) : 0.00918
Pi_(1,0,1) : 0.242
Pi_(1,0,0) : 0.0058
Pi_(0,1,1) : 0.0329
Pi_(0,1,0) : 0.00352
Pi_(0,0,1) : 0.305
Pi_(0,0,0) : 0.0262

EM iteration 32times.