# data analysis and wrangling
import pandas as pd
import numpy as np
import random as rnd
import scipy.stats as stats
from scipy.stats.mstats import winsorize
pd.set_option('display.max_columns', 100)
pd.set_option('display.max_rows', 100)


# visualization
import missingno as msno
import seaborn as sns
import matplotlib.pyplot as plt
%matplotlib inline

# ignore warnings
import warnings
warnings.filterwarnings(action='ignore')


# machine learning
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import LogisticRegression
from sklearn.svm import SVC, LinearSVC
from sklearn.ensemble import RandomForestClassifier
from sklearn.neighbors import KNeighborsClassifier
from sklearn.naive_bayes import GaussianNB
from sklearn.linear_model import Perceptron
from sklearn.linear_model import SGDClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.model_selection import train_test_split

data_2017 = pd.read_csv('./NHIS_2017_2018_100m/NHIS_OPEN_GJ_2017_100.csv')
data_2018 = pd.read_csv('./NHIS_2017_2018_100m/NHIS_OPEN_GJ_2018_100.csv')
df = pd.concat([data_2017,data_2018])

drop_cols = ['HCHK_YEAR', 'IDV_ID', 'HCHK_OE_INSPEC_YN', 'CRS_YN', 'TTH_MSS_YN', 'ODT_TRB_YN', 'TTR_YN', 'WSDM_DIS_YN', 'DATA_STD__DT']
df.drop(drop_cols, axis = 1, inplace = True)

# error시 : sudo conda install missingno
msno.matrix(df)
plt.show()

sns.kdeplot(df.BLDS)
plt.title("BLDS distribution")
plt.show()

p = len(df[df['BLDS'] >= 126])
print('[공복혈당 126 이상]당뇨환자 수치 : {}%'.format(100 * p / len(df)))
p = len(df[(100 <= df['BLDS']) & (df['BLDS'] < 126)])
print('[공복혈당 100 이상 126 미만]전 당뇨환자 수치 : {}%'.format(100 * p / len(df)))
p = len(df[df['BLDS'] < 100])
print('[공복혈당 100 미만]정상인 수치 : {}%'.format(100 * p / len(df)))

df.loc[ df['BLDS'] < 100, 'BLDS'] = 0
df.loc[(100 <= df['BLDS']) & (df['BLDS'] < 126), 'BLDS'] = 1
df.loc[126 <= df['BLDS'], 'BLDS'] = 2
df = df.dropna(subset = ['BLDS'], how = 'any', axis=0)

msno.matrix(df)
plt.show()
print(df.isnull().sum())

mapper = {1 : 'male', 2 : 'female'}
df['SEX'].replace(mapper, inplace=True)
df_ = df[["SEX", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_.groupby(['SEX'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

df['SEX'] = df['SEX'].astype('category')
pd.get_dummies(df['SEX'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

mapper = {11 : 'SEOUL',
          26 : 'BUSAN',
          27 : 'DAGU',
          28 : 'INCHEON',
          29 : 'KWANGJU',
          30 : 'DAJEON',
          31 : 'ULSAN',
          36 : 'SEJONG',
          41 : 'GYEONGGI',
          42 : 'GANGWON',
          43 : 'CB',
          44 : 'CN',
          45 : 'JB',
          46 : 'JN',
          47 : 'GB',
          48 : 'GN',
          49 : 'JEJU',
          50 : 'FOREIGN'}
df['SIDO'].replace(mapper, inplace=True)
df_ = df[["SIDO", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_[["SIDO", "BLDS"]].groupby(['SIDO'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

df['SIDO'] = df['SIDO'].astype('category')
pd.get_dummies(df['SIDO'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

mapper = {1: 'Normal', 2 : 'DEAF'}
df['HEAR_LEFT'].replace(mapper, inplace=True)
df_ = df[["HEAR_LEFT", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_[["HEAR_LEFT", "BLDS"]].groupby(['HEAR_LEFT'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

normal_case = (df['BLDS'] == 1) & (df['HEAR_LEFT'].isnull()) # fill 0
deaf_case = (df['BLDS'] == 0) & (df['HEAR_LEFT'].isnull()) # fill 1
df.loc[normal_case,'HEAR_LEFT'] = df.loc[normal_case,'HEAR_LEFT'].fillna('Normal')
df.loc[deaf_case,'HEAR_LEFT'] = df.loc[deaf_case,'HEAR_LEFT'].fillna('DEAF')

df['HEAR_LEFT'] = df['HEAR_LEFT'].astype('category')
pd.get_dummies(df['HEAR_LEFT'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

mapper = {1: 'Normal', 2 : 'DEAF'}
df['HEAR_RIGHT'].replace(mapper, inplace=True)
df_ = df[["HEAR_RIGHT", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_[["HEAR_RIGHT", "BLDS"]].groupby(['HEAR_RIGHT'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

normal_case = (df['BLDS'] == 1) & (df['HEAR_RIGHT'].isnull()) # fill 0
deaf_case = (df['BLDS'] == 0) & (df['HEAR_RIGHT'].isnull()) # fill 1
df.loc[normal_case,'HEAR_RIGHT'] = df.loc[normal_case,'HEAR_RIGHT'].fillna('Normal')
df.loc[deaf_case,'HEAR_RIGHT'] = df.loc[deaf_case,'HEAR_RIGHT'].fillna('DEAF')

df['HEAR_RIGHT'] = df['HEAR_RIGHT'].astype('category')
pd.get_dummies(df['HEAR_RIGHT'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

mapper = {1: 'No', 2 : 'No', 3: 'Yes'}
df['SMK_STAT_TYPE_CD'].replace(mapper, inplace=True)
df_ = df[["SMK_STAT_TYPE_CD", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_[["SMK_STAT_TYPE_CD", "BLDS"]].groupby(['SMK_STAT_TYPE_CD'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

smoke = (df['BLDS'] == 1) & (df['SMK_STAT_TYPE_CD'].isnull()) # fill 0
non_smoke = (df['BLDS'] == 0) & (df['SMK_STAT_TYPE_CD'].isnull()) # fill 1
df.loc[smoke,'SMK_STAT_TYPE_CD'] = df.loc[smoke,'SMK_STAT_TYPE_CD'].fillna('Yes')
df.loc[non_smoke,'SMK_STAT_TYPE_CD'] = df.loc[non_smoke,'SMK_STAT_TYPE_CD'].fillna('No')

df['SMK_STAT_TYPE_CD'] = df['SMK_STAT_TYPE_CD'].astype('category')
pd.get_dummies(df['SMK_STAT_TYPE_CD'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

mapper = {1: 'Yes', 0 : 'No'}
df['DRK_YN'].replace(mapper, inplace=True)
df_ = df[["DRK_YN", "BLDS"]]
df_['BLDS'] = np.where(df_['BLDS'] == 2, 1, df_['BLDS'])
df_[["DRK_YN", "BLDS"]].groupby(['DRK_YN'], as_index=False).mean().sort_values(by='BLDS', ascending=False)

# 결측치 행 삭제
df = df.dropna(subset = ['DRK_YN'], how = 'any', axis=0)

df['DRK_YN'] = df['DRK_YN'].astype('category')
pd.get_dummies(df['DRK_YN'], drop_first=True)

df = pd.get_dummies(df, drop_first=True)

plt.figure(figsize=(15,10))
sns.countplot(df['AGE_GROUP'])
plt.title("AGE_GROUP",fontsize=15)
plt.show()

plt.figure(figsize=(15,10))
sns.countplot(df['HEIGHT'])
plt.title("HEIGHT",fontsize=15)
plt.show()

plt.figure(figsize=(15,10))
sns.countplot(df['WEIGHT'])
plt.title("WEIGHT",fontsize=15)
plt.show()

df['BMI'] = df['WEIGHT'] / (df['HEIGHT']/100)**2

sns.kdeplot(df.BMI)
plt.title("BMI distribution")
plt.show()

df['PIBW'] = np.where(df['SEX_male'] == 1, (100*df['WEIGHT'])/(22*(df['HEIGHT']/100)**2), (100*df['WEIGHT'])/(21*(df['HEIGHT']/100)**2))

sns.kdeplot(df.PIBW)
plt.title("PIBW distribution")
plt.show()

# 900이상의 허리둘레는 잘못 기입된 값이므로 이상치로 처리한다.
df['WAIST'] = np.where(df['WAIST']>900, np.nan, df['WAIST'])
# 40cm미만의 허리둘레는 잘못 기입된 값이므로 이상치로 처리한다.
df['WAIST'] = np.where(df['WAIST']<40, np.nan, df['WAIST'])

sns.kdeplot(df.WAIST)
plt.title("WAIST distribution")
plt.show()

fig, (ax1,ax2) = plt.subplots(ncols=2)
fig.set_size_inches(12, 5)
sns.scatterplot(data=df, y="BMI", x = 'WAIST', ax=ax1)
sns.scatterplot(data=df, y="PIBW", x = 'WAIST', ax=ax2)
print(df[['PIBW', 'BMI','WAIST']].corr())

df_copy = df[df['WAIST'].notnull()][['BMI', 'PIBW', 'WAIST', 'HEIGHT', 'WEIGHT']]
X_train, X_test, y_train, y_test = train_test_split(df_copy[['BMI', 'PIBW', 'HEIGHT', 'WEIGHT']], 
                                                    df_copy['WAIST'], random_state=42)
lr = LinearRegression()
lr.fit(X_train, y_train)
lr.score(X_test, y_test)

df['WAIST_pred'] = lr.predict(df[['BMI','PIBW', 'HEIGHT', 'WEIGHT']])
df['WAIST'] = np.where(df['WAIST']>0, df['WAIST'], df['WAIST_pred'])

sns.kdeplot(df.WAIST)
plt.title("WAIST distribution")
plt.show()

# 임시 feature WAIST_pred drop
drop_cols = ['WAIST_pred']
df.drop(drop_cols, axis = 1, inplace = True)

# 실명 9.9 -> 0
df['SIGHT_LEFT'] = np.where(df['SIGHT_LEFT']>2.5, 0, df['SIGHT_LEFT'])

plt.figure(figsize=(15,10))
sns.countplot(df['SIGHT_LEFT'])
plt.title("SIGHT_LEFT",fontsize=15)
plt.show()

df['SIGHT_LEFT'].fillna(1, inplace = True)

# 실명 9.9 -> 0
df['SIGHT_RIGHT'] = np.where(df['SIGHT_RIGHT']>2.5, 0, df['SIGHT_RIGHT'])

plt.figure(figsize=(15,10))
sns.countplot(df['SIGHT_RIGHT'])
plt.title("SIGHT_RIGHT",fontsize=15)
plt.show()

df['SIGHT_RIGHT'].fillna(1, inplace = True)

plt.figure(figsize=(15,10))
sns.countplot(df['OLIG_PROTE_CD'])
plt.title("OLIG_PROTE_CD",fontsize=15)
plt.show()

df['OLIG_PROTE_CD'].fillna(1, inplace = True)

df = df.dropna(subset = ['BP_HIGH'], how = 'any', axis=0)

sns.distplot(df.BP_HIGH)
plt.title("BP_HIGH")
plt.show()

df['BP_HIGH_level'] = np.where(160 <= df['BP_HIGH'], 3, df['BP_HIGH'])
df['BP_HIGH_level'] = np.where((140<=df['BP_HIGH']) & (df['BP_HIGH']<160), 2, df['BP_HIGH_level'])
df['BP_HIGH_level'] = np.where((121<=df['BP_HIGH']) & (df['BP_HIGH']<140), 1, df['BP_HIGH_level'])
df['BP_HIGH_level'] = np.where((100<=df['BP_HIGH']) & (df['BP_HIGH']<=120), 0, df['BP_HIGH_level'])
df['BP_HIGH_level'] = np.where((df['BP_HIGH'].min()<=df['BP_HIGH']) & (df['BP_HIGH']<100), -1, df['BP_HIGH_level'])

plt.figure(figsize=(15,10))
sns.countplot(df['BP_HIGH_level'])
plt.title("BP_HIGH_level",fontsize=15)
plt.show()

df = df.dropna(subset = ['BP_LWST'], how = 'any', axis=0)

sns.distplot(df.BP_LWST)
plt.title("BP_LWST")
plt.show()

df['BP_LWST_level'] = np.where(100 <= df['BP_LWST'], 3, df['BP_LWST'])
df['BP_LWST_level'] = np.where((90<=df['BP_LWST']) & (df['BP_LWST']<100), 2, df['BP_LWST_level'])
df['BP_LWST_level'] = np.where((81<=df['BP_LWST']) & (df['BP_LWST']<90), 1, df['BP_LWST_level'])
df['BP_LWST_level'] = np.where((60<=df['BP_LWST']) & (df['BP_LWST']<=80), 0, df['BP_LWST_level'])
df['BP_LWST_level'] = np.where((df['BP_LWST'].min()<=df['BP_LWST']) & (df['BP_LWST']<60), -1, df['BP_LWST_level'])

plt.figure(figsize=(15,10))
sns.countplot(df['BP_LWST_level'])
plt.title("BP_LWST_level",fontsize=15)
plt.show()

df = df.dropna(subset = ['HMG'], how = 'any', axis=0)

sns.distplot(df.HMG)
plt.title("HMG")
plt.show()

sns.boxplot(data=df, x='HMG')

# 4개의 feature 결측치 제거
df = df.dropna(subset = ['CREATININE', 'SGOT_AST', 'SGPT_ALT', 'GAMMA_GTP'], how = 'any', axis=0)

outlier = df[['CREATININE', 'SGOT_AST', 'SGPT_ALT', 'GAMMA_GTP']].copy()
for col in outlier :
  q1 = outlier[col].quantile(0.25)
  q3 = outlier[col].quantile(0.75)
  iqr = q3 - q1
  lower_bound = q1 - (iqr * 1.5)
  upper_bound = q3 + (iqr * 1.5)
  print('{}\'s upper bound : {}, lower bount : {}'.format(col, round(upper_bound, 2), round(lower_bound, 2)))
  outlier[col] = np.where(outlier[col] > upper_bound, np.nan, outlier[col])
  outlier[col] = np.where(outlier[col] < lower_bound, np.nan, outlier[col])
msno.matrix(outlier)
plt.show()

outlier.isnull().sum()
print(outlier.isnull().sum())
print('-'*30)
print('전체 이상치 개수 :', len(outlier) - len(outlier.dropna()))

fig, axes = plt.subplots(nrows = 3, ncols=4, figsize=(18, 12))
for i, col in enumerate(outlier) :
  sns.kdeplot(df[col], ax = axes[0][i])
  sns.countplot(outlier[col],ax = axes[1][i])
  sns.boxplot(data=outlier, x=col, orient="v", ax = axes[2][i])
plt.setp(axes.flat, xlabel=None, ylabel=None)
pad = 5
rows = ['include outlier', 'exclude outlier', 'exclude outlier boxplot']
for ax, row in zip(axes[:,0], rows):
    ax.annotate(row, xy=(0, 0.5), xytext=(-ax.yaxis.labelpad - pad, 0),
                xycoords=ax.yaxis.label, textcoords='offset points',
                size='large', ha='right', va='center')

fig.tight_layout()
plt.show()

# delete outlier
outlier = ['CREATININE', 'SGOT_AST', 'SGPT_ALT', 'GAMMA_GTP']
for col in outlier :
  q1 = df[col].quantile(0.25)
  q3 = df[col].quantile(0.75)
  iqr = q3 - q1
  lower_bound = q1 - (iqr * 1.5)
  upper_bound = q3 + (iqr * 1.5)
  df[col] = np.where(df[col] > upper_bound, np.nan, df[col])
  df[col] = np.where(df[col] < lower_bound, np.nan, df[col])

df = df.dropna(subset = outlier, how = 'any', axis=0)

msno.matrix(df)
plt.show()
print(df.isnull().sum())

outlier = df[['TOT_CHOLE', 'TRIGLYCERIDE', 'HDL_CHOLE', 'LDL_CHOLE']].copy()
a = len(outlier.dropna())
for col in outlier :
  q1 = outlier[col].quantile(0.25)
  q3 = outlier[col].quantile(0.75)
  iqr = q3 - q1
  lower_bound = q1 - (iqr * 1.5)
  upper_bound = q3 + (iqr * 1.5)
  print('{}\'s upper bound : {}, lower bount : {}'.format(col, round(upper_bound, 2), round(lower_bound, 2)))
  outlier[col] = np.where(outlier[col] > upper_bound, np.nan, outlier[col])
  outlier[col] = np.where(outlier[col] < lower_bound, np.nan, outlier[col])
msno.matrix(outlier)
plt.show()

outlier.isnull().sum()
print(outlier.isnull().sum())
print('-'*30)
print('전체 이상치 및 outlier 개수 :', len(outlier) - len(outlier.dropna()))

fig, axes = plt.subplots(nrows = 3, ncols=4, figsize=(18, 12))
for i, col in enumerate(outlier) :
  sns.kdeplot(df[col], ax = axes[0][i])
  sns.countplot(outlier[col],ax = axes[1][i])
  sns.boxplot(data=outlier, x=col, orient="v", ax = axes[2][i])
plt.setp(axes.flat, xlabel=None, ylabel=None)

pad = 5
rows = ['include outlier', 'exclude outlier', 'exclude outlier boxplot']
for ax, row in zip(axes[:,0], rows):
    ax.annotate(row, xy=(0, 0.5), xytext=(-ax.yaxis.labelpad - pad, 0),
                xycoords=ax.yaxis.label, textcoords='offset points',
                size='large', ha='right', va='center')

fig.tight_layout()
plt.show()

# 이상치 제거
cols = ['TOT_CHOLE', 'TRIGLYCERIDE', 'HDL_CHOLE', 'LDL_CHOLE']
for col in cols :
  q1 = df[col].quantile(0.25)
  q3 = df[col].quantile(0.75)
  iqr = q3 - q1
  lower_bound = q1 - (iqr * 1.5)
  upper_bound = q3 + (iqr * 1.5)
  df = df[(lower_bound<=df[col]) & (df[col] <= upper_bound) | (df[col].isnull())]

countinous = ['AGE_GROUP', 'HEIGHT', 'WEIGHT', 'WAIST', 'SIGHT_LEFT', 
              'SIGHT_RIGHT', 'BP_HIGH', 'BP_LWST', 'BLDS', 'TOT_CHOLE',
              'TRIGLYCERIDE', 'HDL_CHOLE', 'LDL_CHOLE', 'HMG', 'OLIG_PROTE_CD',
              'CREATININE', 'SGOT_AST', 'SGPT_ALT', 'GAMMA_GTP', 'BMI', 'PIBW', 
              'BP_HIGH_level', 'BP_LWST_level']
corrMatt = df[countinous].corr()
mask = np.array(corrMatt)
mask[np.tril_indices_from(mask)] = False
f,ax = plt.subplots(figsize=(30, 30))
#sns.heatmap(corrMatt, annot=True, linewidths=.5, fmt= '.1f',ax=ax)
sns.heatmap(corrMatt, mask=mask,vmax=.8, square=True,annot=True, ax = ax)

from sklearn.linear_model import LinearRegression, Ridge, Lasso, ElasticNet
from sklearn.preprocessing import StandardScaler


data = df.copy()
data.dropna(inplace = True)
y_tot = data['TOT_CHOLE']
y_tri = data['TRIGLYCERIDE']
y_hdl = data['HDL_CHOLE']
y_ldl = data['LDL_CHOLE']

drop_cols = ['TOT_CHOLE', 'TRIGLYCERIDE', 'HDL_CHOLE', 'LDL_CHOLE', 'BLDS']
data.drop(drop_cols, axis = 1, inplace = True)

X_train, X_test, y_tot_train, y_tot_test = train_test_split(data, y_tot, test_size = 0.3, random_state=42)  
X_train, X_test, y_tri_train, y_tri_test = train_test_split(data, y_tri, test_size = 0.3, random_state=42)  
X_train, X_test, y_hdl_train, y_hdl_test = train_test_split(data, y_hdl, test_size = 0.3, random_state=42)  
X_train, X_test, y_ldl_train, y_ldl_test = train_test_split(data, y_ldl, test_size = 0.3, random_state=42)  

sc = StandardScaler().fit(X_train)
X_train_sc = sc.transform(X_train)
X_test_sc = sc.transform(X_test)

model_tot = LinearRegression()
model_tot.fit(X_train_sc, y_tot_train)
print("TOT score : ",model_tot.score(X_test_sc, y_tot_test))
y_pred = model_tot.predict(X_test_sc) 
print("TOT MSE : ", sum((y_tot_test - y_pred)**2) / len(y_tot_test))

model_tri = LinearRegression()
model_tri.fit(X_train_sc, y_tri_train)
print("TRI score : ",model_tri.score(X_test_sc, y_tri_test))
y_pred = model_tri.predict(X_test_sc) 
print("TRI MSE : ", sum((y_tri_test - y_pred)**2) / len(y_tri_test))

model_hdl = LinearRegression()
model_hdl.fit(X_train_sc, y_hdl_train)
print("HDL score : ",model_hdl.score(X_test_sc, y_hdl_test))
y_pred = model_hdl.predict(X_test_sc) 
print("HDL MSE : ", sum((y_hdl_test - y_pred)**2) / len(y_hdl_test))

model_ldl = LinearRegression()
model_ldl.fit(X_train_sc, y_ldl_train)
print("LDL score : ",model_ldl.score(X_test_sc, y_ldl_test))
y_pred = model_ldl.predict(X_test_sc) 
print("LDL MSE : ", sum((y_ldl_test - y_pred)**2) / len(y_ldl_test))

df.dropna(inplace = True)
msno.matrix(df)
plt.show()

df.describe()

df.shape

fat = (df['WAIST'] >= 90) & (df['SEX_male'] == 1) | (df['WAIST'] >= 85) & (df['SEX_male'] == 0)
df['ABD_FAT'] = np.where(fat, 1, 0)

warn = (df['GAMMA_GTP'] > 63) & (df['SEX_male'] == 1) | (df['GAMMA_GTP'] > 35) & (df['SEX_male'] == 0)
df['GAMMA_GTP_level'] = np.where(warn, 1, 0)