[심화세션]_4주차

주교재

[7장] 물류 네트워크 최적 설계를 위한 테크닉 10

61. 운송 최적화 문제를 풀어보자

# 초기 설정  #
np.random.seed(1)
nw = len(df_tc.index)
nf = len(df_tc.columns)
pr = list(product(range(nw), range(nf)))

# 수리 모델 작성 #

# 모델 생성
m1 = model_min()
# v1에 dict 형식으로 LpVariable 정의 -> 'v1_2' 이런 형식으로
v1 = {(i,j):LpVariable('v%d_%d'%(i,j),lowBound=0) for i,j in pr}

# 목적함수 정의 : (거래 비용 * 거래량)의 합
m1 += lpSum(df_tc.iloc[i][j]*v1[i,j] for i,j in pr)
for i in range(nw):
    m1 += lpSum(v1[i,j] for j in range(nf)) <= df_supply.iloc[0][i]
for j in range(nf):
    m1 += lpSum(v1[i,j] for i in range(nw)) >= df_demand.iloc[0][j]
# 최적화 문제 해결
m1.solve()

# 총 운송 비용 계산 #
df_tr_sol = df_tc.copy() # df_tr_sol에 df_tc 복사해서 저장
total_cost = 0
for k,x in v1.items():
    i,j = k[0],k[1]
    df_tr_sol.iloc[i][j] = value(x)
    total_cost += df_tc.iloc[i][j]*value(x)
    
print(df_tr_sol)
print("총 운송 비용:"+str(total_cost))

• LpVariable: Pulp, 결정변수 정의, 최적화 문제서 각 변수의 값을 어떻게 최적화할지 모델링
  • 주요 매개변수: name, lowBound, upBound, cat(카테고리, ‘Continuous’ or ‘Integer’)
• value(): Pulp, 선형 계획 문제서 최적화된 변수의 값 추출, solve로 해결한 뒤에 호출
• df.iloc[row_index, column_index]: pandas, 위치 기반 인덱싱, 필터링 불가, 순차적 접근 시에 유용
• df.loc[row_label, column_label]: pandas, 레이블 기반 인덱싱, 필터링 가능, 이름으로 접근해야할 때 유용

62. 최적 운송 경로를 네트워크로 확인하자

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx

# 데이터 불러오기
df_tr = df_tr_sol.copy()
df_pos = pd.read_csv('trans_route_pos.csv')

# 객체 생성
G = nx.Graph()

# 노드 설정: df의 모든 칼럼에 대해 노드 설정
for i in range(len(df_pos.columns)):
    G.add_node(df_pos.columns[i])

# 엣지 설정 & 엣지의 가중치 리스트화
num_pre = 0
edge_weights = []
size = 0.1 # 가중치 스케일 조절(엣지 두께 조절)
for i in range(len(df_pos.columns)):
    for j in range(len(df_pos.columns)):
        # 자기 자신으로 이어지는 것이 아니면 엣지 추가함: 가능한 모든 경우의 엣지 생성
        if not (i==j):
            # 엣지 추가
            G.add_edge(df_pos.columns[i],df_pos.columns[j])
            # 엣지 가중치 추가
            if num_pre<len(G.edges):
                num_pre = len(G.edges)
                weight = 0

                # if문이 True면(0만 아니면), 이어지는 코드 블럭 실행, df가 대칭이 아닐 경우에 대비해 elif문에서 i, j 위치를 바꾸어 한번 더 실행 
                if (df_pos.columns[i] in df_tr.columns)and(df_pos.columns[j] in df_tr.index):
                    if df_tr[df_pos.columns[i]][df_pos.columns[j]]:
                        weight = df_tr[df_pos.columns[i]][df_pos.columns[j]]*size
                elif(df_pos.columns[j] in df_tr.columns)and(df_pos.columns[i] in df_tr.index):
                    if df_tr[df_pos.columns[j]][df_pos.columns[i]]:
                        weight = df_tr[df_pos.columns[j]][df_pos.columns[i]]*size
                edge_weights.append(weight)
                

# 좌표 설정
pos = {}
for i in range(len(df_pos.columns)):
    node = df_pos.columns[i]
    pos[node] = (df_pos[node][0],df_pos[node][1])
    
# 그리기
nx.draw(G, pos, with_labels=True,font_size=16, node_size = 1000, node_color='k', font_color='w', width=edge_weights)

# 표시
plt.show()

• G.add_node(n): 노드 n을 그래프에 추가
• G.add_edge(u,v): u와 v 사이에 엣지를 추가

63. 최적운송경로가 제약조건을 만족하는지 확인하자

# 제약조건 계산함수
# 수요측
def condition_demand(df_tr,df_demand):
    # 모든 요소 0으로 설정
    flag = np.zeros(len(df_demand.columns))
    for i in range(len(df_demand.columns)):
        temp_sum = sum(df_tr[df_demand.columns[i]])
        # df_tr의 수요합이 df_demand의 수요량보다 크면 flag = 1 설정
        if (temp_sum>=df_demand.iloc[0][i]):
            flag[i] = 1
    return flag
            
# 공급측
def condition_supply(df_tr,df_supply):
    # 모든 요소 0으로 설정
    flag = np.zeros(len(df_supply.columns))
    for i in range(len(df_supply.columns)):
        temp_sum = sum(df_tr.loc[df_supply.columns[i]])
        # df_tr의 공급합이 df_supply의 공급량보다 적으면 flag = 1 설정
        if temp_sum<=df_supply.iloc[0][i]:
            flag[i] = 1
    return flag

print("수요 조건 계산 결과:"+str(condition_demand(df_tr_sol,df_demand)))
print("공급 조건 계산 결과:"+str(condition_supply(df_tr_sol,df_supply)))

# 수요 조건 계산 결과:[1. 1. 1. 1.]
# 공급 조건 계산 결과:[1. 1. 1.]

64. 생산계획 데이터를 불러오자

지금까지는 운송 비용 최적화를 계산했다. 이제부턴 생산 계획 최적화를 계산해보자.

65. 이익을 계산하는 함수를 만들자

생산계획 최적화

1. 목적함수, 제약 조건 정의
2. 제약 조건 아래서 목적함수 최소(대)화 변수 조합 찾기

# 이익 계산 함수
def product_plan(df_profit,df_plan):
    profit = 0
    for i in range(len(df_profit.index)):
        for j in range(len(df_plan.columns)):
            profit += df_profit.iloc[i][j]*df_plan.iloc[i][j]
    return profit

print("총 이익:"+str(product_plan(df_profit,df_plan)))

66. 생산최적화 문제를 풀어보자

import pandas as pd
from pulp import LpVariable, lpSum, value
from ortoolpy import model_max, addvars, addvals


df = df_material.copy()
inv = df_stock

m = model_max()
v1 = {(i):LpVariable('v%d'%(i),lowBound=0) for i in range(len(df_profit))}
m += lpSum(df_profit.iloc[i]*v1[i] for i in range(len(df_profit)))
for i in range(len(df_material.columns)):
    m += lpSum(df_material.iloc[j,i]*v1[j] for j in range(len(df_profit)) ) <= df_stock.iloc[:,i]
m.solve()

df_plan_sol = df_plan.copy()
for k,x in v1.items():
    df_plan_sol.iloc[k] = value(x)
print(df_plan_sol)
print("총 이익:"+str(value(m.objective)))

67. 최적생산계획이 제약조건을 만족하는지 확인하자

최적화 문제를 풀 때에 최적화 계산을 한 결과를 이해하지 않고 그대로 받아들여서는 안된다!

# 제약 조건 계산 함수
def condition_stock(df_plan,df_material,df_stock):
    flag = np.zeros(len(df_material.columns))
    for i in range(len(df_material.columns)):  
        temp_sum = 0
        for j in range(len(df_material.index)):  
            temp_sum = temp_sum + df_material.iloc[j][i]*float(df_plan.iloc[j])
        if (temp_sum<=float(df_stock.iloc[0][i])):
            flag[i] = 1
        print(df_material.columns[i]+"  사용량:"+str(temp_sum)+", 재고:"+str(float(df_stock.iloc[0][i])))
    return flag

print("제약 조건 계산 결과:"+str(condition_stock(df_plan_sol,df_material,df_stock)))

68. 물류 네트워크 설계 문제를 풀어보자

지금까지는 운송 경로와 생산 계획 최적화 문제를 따로 생각했으나 실제로는 이 둘을 동시에 고려해야 한다.

목적함수를 ‘운송 비용+ 제조 비용’으로 정의하고,

제약조건은 각 대리점의 판매 수가 수요 수를 넘는 것으로 정의하자.

import numpy as np
import pandas as pd

# 각 변수들에 해당하는 리스트, 튜플 생성
제품 = list('AB')
대리점 = list('PQ')
공장 = list('XY')
레인 = (2,2)

# 운송비 #
tbdi = pd.DataFrame(((j,k) for j in 대리점 for k in 공장), columns=['대리점','공장'])
tbdi['운송비'] = [1,2,3,1]
print(tbdi)

# 수요 #
tbde = pd.DataFrame(((j,i) for j in 대리점 for i in 제품), columns=['대리점','제품'])
tbde['수요'] = [10,10,20,20]
print(tbde)

# 생산 #
tbfa = pd.DataFrame(((k,l,i,0,np.inf) for k,nl in zip (공장,레인) for l in range(nl) for i in 제품), 
                    columns=['공장','레인','제품','하한','상한'])
tbfa['생산비'] = [1,np.nan,np.nan,1,3,np.nan,5,3]
tbfa.dropna(inplace=True)
tbfa.loc[4,'상한']=10
print(tbfa)

from ortoolpy import logistics_network
# tbde, tbdi, tbfa를 입력으로 사용, 사용할 칼럼 이름 지정
# logistics_network()는 여러 개의 값을 반환하는데, 이 중에서 tbdi2만 관심 있으므로 나머지는 무시하기 위해, 두 번째 값만 사용하기 위해 ' _, tbdi2, _'로 표현하고 이름 짓기
_, tbdi2, _ = logistics_network(tbde, tbdi, tbfa,dep = "대리점", dem = "수요",fac = "공장",prd = "제품",tcs = "운송비",pcs = "생산비",lwb = "하한",upb = "상한")

print(tbfa)
print(tbdi2)

• zip(공장, 레인): 파이썬, 여러 리스트를 병렬로 묶어주기 위한 함수

예시:

list1 = [1, 2, 3]
list2 = ['a', 'b', 'c']

zipped = zip(list1, list2)
print(list(zipped))

69. 최적 네트워크의 운송비용과 그 내역을 계산하자

# 새로운 df형태로 지정
tbdi2 = tbdi2[["공장","대리점","운송비","제품","VarX","ValX"]]
tbdi2

trans_cost = 0
for i in range(len(tbdi2.index)):
    trans_cost += tbdi2["운송비"].iloc[i]*tbdi2["ValX"].iloc[i]
print("총 운송비:"+str(trans_cost))

70. 최적네트워크의 생산비용과 그 내역을 계산하자

tbfa

product_cost = 0
for i in range(len(tbfa.index)):
    product_cost += tbfa["생산비"].iloc[i]*tbfa["ValY"].iloc[i]
print("총 생산비:"+str(product_cost))

[8장] 수치 시뮬레이션으로 소비자의 행동을 예측하는 테크닉 10

71. 인간관계 네트워크를 가시화해보자

import networkx as nx
import matplotlib.pyplot as plt

# 그래프 객체 생성
G = nx.Graph()

# 노드 설정
NUM = len(df_links.index)
for i in range(1,NUM+1):
    node_no = df_links.columns[i].strip("Node")
    G.add_node(str(node_no))

# 엣지 설정
for i in range(NUM):
    for j in range(NUM):
        node_name = "Node" + str(j)
        if df_links[node_name].iloc[i]==1:
            G.add_edge(str(i),str(j))
        
# 그리기
nx.draw_networkx(G,node_color="k", edge_color="k", font_color="w")
plt.show()

• nx.draw_networkx(): nx.draw()보다 복잡한 그래프 시각화 가능, 자동으로 연결이 많은 노드를 중심에 오도록 설정

72. 입소문에 의한 정보 전파 모습을 가시화해보자

# 랜덤으로 연결여부 결정
def determine_link(percent):
    rand_val = np.random.rand()
    if rand_val<=percent:
        return 1
    else:
        return 0

# 노드가 active 상태이면, j와의 연결 활성화 결정
def simulate_percolation(num, list_active, percent_percolation):
    for i in range(num):
        if list_active[i]==1:
            for j in range(num):
                node_name = "Node" + str(j)
                if df_links[node_name].iloc[i]==1:
                    if determine_link(percent_percolation)==1:
                        list_active[j] = 1
    return list_active

percent_percolation = 0.1
T_NUM = 36
NUM = len(df_links.index)
list_active = np.zeros(NUM)
list_active[0] = 1

list_timeSeries = []
for t in range(T_NUM):
    list_active = simulate_percolation(NUM, list_active, percent_percolation)
    # list_active가 매 시간 단계에서 변경되기에 이전 상태를 보존하기 위해 .copy() 사용
    list_timeSeries.append(list_active.copy())

# 액티브 노드 가시화 #
def active_node_coloring(list_active):
    #print(list_timeSeries[t])
    list_color = []
    for i in range(len(list_timeSeries[t])):
        if list_timeSeries[t][i]==1:
            list_color.append("r")
        else:
            list_color.append("k")
    #print(len(list_color))
    return list_color

# 그리기
t = 0
nx.draw_networkx(G,font_color="w",node_color=active_node_coloring(list_timeSeries[t]))
plt.show()

t=11,35로 가시화해보면, 12개월까지는 완만한 전파였지만, 오랜 시간이 흐른 뒤 거의 전원에게 전파되는 모습을 볼 수 있다.

73. 입소문 수의 시계열 변화를 그래프화 해보자

# 시계열 그래프 그리기
list_timeSeries_num = []
# list_timeSeries는 percolation함수에서 노드 활성 시간을 저장한 리스트
for i in range(len(list_timeSeries)):
    list_timeSeries_num.append(sum(list_timeSeries[i]))

plt.plot(list_timeSeries_num)
plt.show()

74. 회원수의 시계열 변화를 시뮬레이션해보자

def simulate_population(num, list_active, percent_percolation, percent_disapparence,df_links):
    # 확산 #
    for i in range(num):
        if list_active[i]==1:
            for j in range(num):
                if df_links.iloc[i][j]==1:
                    if determine_link(percent_percolation)==1:
                        list_active[j] = 1
    # 소멸 #
    for i in range(num):
        if determine_link(percent_disapparence)==1:
            list_active[i] = 0
    return list_active

percent_percolation = 0.1
percent_disapparence = 0.05
T_NUM = 100
NUM = len(df_links.index)
list_active = np.zeros(NUM)
list_active[0] = 1

list_timeSeries = []
for t in range(T_NUM):
    list_active = simulate_population(NUM, list_active, percent_percolation, percent_disapparence,df_links)
    list_timeSeries.append(list_active.copy())

# 시계열 그래프 그리기
list_timeSeries_num = []
for i in range(len(list_timeSeries)):
    list_timeSeries_num.append(sum(list_timeSeries[i]))

plt.plot(list_timeSeries_num)
plt.show()

percent_disapparence = 0.2
list_active = np.zeros(NUM)
list_active[0] = 1
list_timeSeries = []
for t in range(T_NUM):
    list_active = simulate_population(NUM, list_active, percent_percolation, percent_disapparence,df_links)
    list_timeSeries.append(list_active.copy())

# 시계열 그래프 그리기
list_timeSeries_num = []
for i in range(len(list_timeSeries)):
    list_timeSeries_num.append(sum(list_timeSeries[i]))

plt.plot(list_timeSeries_num)
plt.show()

75. 파라미터 전체를 ‘상관관계’를 보면서 파악해보자

# 상관관계 계산
print("상관관계 계산시작")
T_NUM = 100
NUM_PhaseDiagram = 20
phaseDiagram = np.zeros((NUM_PhaseDiagram,NUM_PhaseDiagram))
for i_p in range(NUM_PhaseDiagram):
    for i_d in range(NUM_PhaseDiagram):
        percent_percolation = 0.05*i_p
        percent_disapparence = 0.05*i_d
        list_active = np.zeros(NUM)
        list_active[0] = 1
        for t in range(T_NUM):
            list_active = simulate_population(NUM, list_active, percent_percolation, percent_disapparence,df_links)
        phaseDiagram[i_p][i_d] = sum(list_active)
print(phaseDiagram)

# 표시
# phaseDiagram 배열을 히트맵으로 표시
plt.matshow(phaseDiagram)
# 색상 막대 속성 표시, 크기 조절
plt.colorbar(shrink=0.8)
plt.xlabel('percent_disapparence')
plt.ylabel('percent_percolation')
# 눈금 위치, 레이블 설정
plt.xticks(np.arange(0.0, 20.0,5), np.arange(0.0, 1.0, 0.25))
plt.yticks(np.arange(0.0, 20.0,5), np.arange(0.0, 1.0, 0.25))
# 축의 눈금 표시 비활성화
plt.tick_params(bottom=False,
                left=False,
                right=False,
                top=False)
plt.show()

76. 실제 데이터를 불러와보자

import pandas as pd

df_mem_links = pd.read_csv("links_members.csv")
df_mem_info = pd.read_csv("info_members.csv")
df_mem_links

77. 링크 수의 분포를 가시화해보자

네트워크의 구조

• 스몰 월드형: 몇 안 되는 스텝으로 전원이 연결됨
• 스케일 프리형: 소수의 연결을 많이 가지는 사람이 허브가 되

NUM = len(df_mem_links.index)
array_linkNum = np.zeros(NUM)
for i in range(NUM):
    array_linkNum[i] = sum(df_mem_links["Node"+str(i)])

plt.hist(array_linkNum, bins=10,range=(0,250))
plt.show()

정규분포와 유사한 형태의 히스토그램을 나타내는 것을 볼 수 있다.

스케일프리형이라면 이 분포가 ‘멱 법칙’에 가까워지며, 링크를 많이 가진 허브가 작동하지 않으면(입소문을 퍼뜨리지 않으면) 입소문이 중간에 퍼지지 않는 특징을 가지고 있다.

그러나 이 분포는 거의 모든 노드가 어느 정도의 링크 수를 가지고 있기에 ‘급격히 입소문이 퍼지지 않는’ 대신에 ‘허브에 의존하지 않고 입소문이 퍼지기 쉽다’라고 할 수 있다.

78. 시뮬레이션을 위해 실제 데이터로부터 파라미터를 추정하자

NUM = len(df_mem_info.index)
# 첫 번째 열은 노드 이름이기에 -1 수행
T_NUM = len(df_mem_info.columns)-1
# 소멸 확률 추정 #
count_active = 0
count_active_to_inactive = 0
for t in range(1,T_NUM):
    for i in range(NUM):
        if (df_mem_info.iloc[i][t]==1):
            count_active_to_inactive += 1
            # 첫 번째 열을 빼고 시작해서 t+1로 시행
            if (df_mem_info.iloc[i][t+1]==0):
                count_active += 1
estimated_percent_disapparence = count_active/count_active_to_inactive

# 확산 확률 추정: 어렵다.. #
count_link = 0
count_link_to_active = 0
count_link_temp = 0
for t in range(T_NUM-1):
    df_link_t = df_mem_info[df_mem_info[str(t)]==1]
    temp_flag_count = np.zeros(NUM)
    for i in range(len(df_link_t.index)):
        df_link_temp = df_mem_links[df_mem_links["Node"+str(df_link_t.index[i])]==1]
        for j in range(len(df_link_temp.index)):
            if (df_mem_info.iloc[df_link_temp.index[j]][t]==0):
                if (temp_flag_count[df_link_temp.index[j]]==0):
                    count_link += 1
                if (df_mem_info.iloc[df_link_temp.index[j]][t+1]==1):
                    if (temp_flag_count[df_link_temp.index[j]]==0):
                        temp_flag_count[df_link_temp.index[j]] = 1 
                        count_link_to_active += 1
estimated_percent_percolation = count_link_to_active/count_link

estimated_percent_disapparence
# 0.10147163541419416

estimated_percent_percolation
# 0.025184661323275185

79. 실제 데이터와 시뮬레이션을 비교하자

percent_percolation = 0.025184661323275185
percent_disapparence = 0.10147163541419416
T_NUM = 24
NUM = len(df_mem_links.index)
list_active = np.zeros(NUM)
list_active[0] = 1
list_timeSeries = []
for t in range(T_NUM):
    list_active = simulate_population(NUM, list_active, percent_percolation, percent_disapparence,df_mem_links)
    list_timeSeries.append(list_active.copy())

list_timeSeries_num = []
for i in range(len(list_timeSeries)):
    list_timeSeries_num.append(sum(list_timeSeries[i]))

T_NUM = len(df_mem_info.columns)-1
list_timeSeries_num_real = []
for t in range(0,T_NUM):
    list_timeSeries_num_real.append(len(df_mem_info[df_mem_info[str(t)]==1].index))

plt.plot(list_timeSeries_num, label = 'simulated')
plt.plot(list_timeSeries_num_real, label = 'real')
plt.xlabel('month')
plt.ylabel('population')
plt.legend(loc='lower right')
plt.show()

80. 시뮬레이션으로 미래를 예측해보자

percent_percolation = 0.025184661323275185
percent_disapparence = 0.10147163541419416
T_NUM = 36
NUM = len(df_mem_links.index)
list_active = np.zeros(NUM)
list_active[0] = 1
list_timeSeries = []
for t in range(T_NUM):
    list_active = simulate_population(NUM, list_active, percent_percolation, percent_disapparence,df_mem_links)
    list_timeSeries.append(list_active.copy())

list_timeSeries_num = []
for i in range(len(list_timeSeries)):
    list_timeSeries_num.append(sum(list_timeSeries[i]))

plt.plot(list_timeSeries_num, label = 'simulated')
plt.xlabel('month')
plt.ylabel('population')
plt.legend(loc='lower right')
plt.show()

부교재

[8장] 통계 모형화

8.1 선형회귀 원리의 확장

선형회귀의 확장 방향성

1. 설명변수의 개수 늘리기 / 유형 변경하기
2. 반응변수의 유형 변경하기
3. 회귀 모형 형태 변경하기

• 다중회귀: 설명변수가 여러 개인 것

• 편회귀계수(partial regression coefficient): 다중회귀에서의 기울기 (b_1, b_2 등)

• 회귀평면: 설명변수가 2개일 때 그린 그래프

cf. 설명변수가 3개 이상이면 그래프로 나타내기 어렵다

회귀분석에서 구한 편회귀계수는 설명변수의 데이터 퍼짐 정도나 단위에 따라 크게 달라지기 때문에 편회귀계수끼리 비교할 수는 없다.

여기서 편회귀계수를 비교하기 위해 표준화편회귀계수를 이용한다. 표준화회귀계수는 회귀분석을 시행하기 전에 각각의 설명변수를 평균 0, 표준편차 1로 변환한 다음, 회귀분석을 시행하여 구한 회귀계수다.

이 표준화편회귀계수는 각 설명변수가 표준편차 단위에서 1 늘었을 때 반응변수의 증감을 나타낸다.

범주형 변수를 설명변수: 범주에는 대소 관계가 없으므로, 회귀분석의 설명변수로 이용할 때는 0 또는 1과 같은 가변수(dummy variable)를 설명변수로 이용하는 등의 요령이 필요하다.

공분산분석(ANCOVA, analysis of covariance): 일반적인 분산분석에 사용하는 데이터와 함께 양적 변수 데이터가 있는 경우에 후보가 되는 방법. 이때 새로 추가한 양적 변수를 공변량(covariate)라고 한다.

• 공분산분석의 사용 조건
  1. 집단 간 회귀의 기울기가 서로 다르지 않을 것 (회귀직선이 평행)
  2. 회귀계수가 0이 아니어야 한다.

고차원 데이터 문제

1. 차원이 늘어날수록 파라미터 추정에 필요한 데이터 양이 폭발적으로 증가, 차원의 저주
2. 다중공선성 문제가 일어나기 쉬워, 모형의 추정 정밀도가 떨어진다.

• 다중공선성(multicollinearity): 설명변수가 여러 개인 다중회귀에서 설명변수 사이에 강한 상관 존재

• 분산팽창인수 VIF(variance inflation factor): 각 설명변수마다 산출되는 다중공선성 정도 측정 지표

• VIF>10: 2개 사이의 상관이 아주 강한 것. 다중공선성 정도 판단에 대한 하나의 기준

다중공선성이 강하다고 판단했을 때는, 서로 상관이 있는 2개 변수 중 하나를 없애거나, 주성분분석 등의 차원 축소 방법을 이용하여 설명변수의 개수를 줄이는 것이 좋다.

8.2 회귀모형의 형태 바꾸기

• 상호작용: 설명변수 간의 상승효과, 선형회귀모형 안에서 c x_i x_j로 도입할 수 있다.

• 상호작용항 추가의 단점
  • 상호작용을 넣으면 해석이 어려워짐
  • 설명변수의 개수가 늘면 상호작용항의 수가 폭발적으로 증가
  • 상호작용의 형태는 다양한데도 곱셈으로만 나타남
  • 설명변수와 상호작용항의 다중공선성 문제
• 상호작용항 적용이 적절한 경우
  • 상호작용이 있다는 것이 선형 연구에서 밝혀졌거나 기대되는 때
  • 데이터에 분명한 상호작용이 있을 때
  • 상호작용 유무에 관심이 있을 때

• 다원배치 분산분석: 여러 개의 요인을 동시에 고려
  • 이원배치 분산분석: 2개의 요인이 있는 경우

비선형회귀: 일반적으로는 선형모델을 사용하므로, 왜 파라미터가 비선형인 모형을 사용하는지에 대한 합리적인 이유가 있어야만 한다.

8.3 일반화선형모형의 개념

• 통계 모형화: 데이터 성질을 고려하면서 확률 모형을 가정하고, 파라미터를 추정하여 모형을 평가하는 일련의 작업

값이 2개인 반응변수 데이터나 음이 아닌 정수인 반응변수로 구성된 데이터에는, 잔차를 이용해서 거리로 모형과 데이터의 적합도를 측정하기보다는, 데이터가 확률적으로 생성되었다고 간주, ‘확률적으로 얼마나 나타나기 쉬운가’에 기반해 데이터에 잘 들어맞는지 평가하는 편이 좋을 수 있다.

여기서 좌변항을 데이터 x에 대한 파라미터 θ의 함수로 본 것을 가능도라고 한다. 가능도가 큰 것은 그 θ에서 얻은 데이터가 나타나기 쉽다는 것을 뜻한다.

• 최대가능도 방법/추정: 가능도를 최대화하는 θ를 찾아서, 데이터에 가장 잘 들어맞는 파라미터 θ를 정하는 방법

로지스틱 회귀: 범주 하나가 일어날 확률을 p로 두고, 설명변수 x가 바뀌었을 때, p가 얼마나 달라지는지를 조사

오즈비

포아송 회귀

다양한 일반화선형모형

• 일반화선형혼합모형

8.4 통계 모형의 평가와 비교

왈드 검정

• 왈드 통계량: 최대가능도 방법으로 얻은 추정값/표준오차

가능도비 검정

AIC (아카이케 정보기준, Akaike information criterion)

• AIC는 새롭게 얻을 데이터를 얼마나 잘 예측할 수 있는지를 바탕으로 모형을 좋음(적합도)을 결정하는 지표다.
• AIC는 실제 데이터에 잘 들어맞는지는 물론, 모형의 파라미터 개수까지 고려하여 모형을 평가함으로써 과대적합이 없는 모형을 선택하는 지표다
• 주의사항은 새롭게 얻을 데이터의 예측도를 높이는 모형을 고르는 것이 목적인 지표이므로, AIC를 최소화한다고 해서 그것이 반드시 실제 모형이지는 않을 수도 있다는 점이다.

BIC (베이즈 정보기준, Bayesian information criterion)

• 최소화하는 것이 좋은 모형
• 표본크기 n에 따라 달라짐, n을 무한대로 하면 실제 모형을 고를 확률이 1이 된다는 일치성 존재

그 밖의 정보기준

• AIC, BIC 이외에도 AICc, DIC, WAIC, WBIC 등을 사용할 수도 있다.

[9장] 가설검정의 주의점

9.1 재현성

• 재현성(reproducibility, replication): 누가 언제 어디서 실험하더라도, 조건이 동일하다면 동일한 결과를 얻을 수 있어야 한다는 것.

• 재현성 위기: 다른 연구자가 동일한 방법과 조건으로 추시했을 때, 같은 결과를 얻지 못했다는 것, 재현성이 없다는 것은 원래 논문의 주장이 잘못되었을 가능성이 있다는 것.

재현 불가능한 원인은?

1. 실험 조건을 동일하게 조성하기 어렵다.
2. 가설검정의 사용 방법
   • p-hacking: p값이 0.05보다 작아지게 조작하는 것

9.2 가설검정의 문제점

• p값의 정의: 귀무가설이 옳다고 가정할 때 실제 관찰한 데이터 이상으로 극단적인 값을 얻을 확률

이 값이 작으면, 귀무가설과 관찰한 데이터 사이에 괴리가 크다는 것을 뜻하며, 아예 유의수준 α를 밑도는 때에는 귀무가설을 기각하는 판단을 내리게 된다.

피셔류 검정과 네이만-피어슨류 검정

• 피셔류 검정: p값을 계산하고, 귀무가설과 관찰한 값의 괴리 정도를 평가한다.
• 네이만-피어슨류 검정: p값이 유의수준 α미만인가 이상인가에만 주목하여, 가설 기각/채택이라는 결과를 내린다.

효과크기 (effect size)

• Cohen’s d: 모집단의 분산을 기준으로 2개의 모집단평균이 얼마나 떨어져 있는지를 나타냄

• 메타 분석: 어떤 현상을 보고한 여러 논문을 통합하여, 결과를 종합적으로 평가하는 방법

• 다양한 효과크기

베이즈 인수 (Bayes factor)

논문이 옳지 않을 확률

• 허위발견율(FDR, false discovery rate): 옳다고 주장된 것 중 위양성인 것의 비율

좋은 가설 세우기

9.3 p-해킹

• p-hacking: 의도와 무관하게, p값을 원하는 방향으로 조작하는 행위
  • 결과를 보며 표본크기를 늘려서는 안됨
  • 마음에 드는 해석만 보고해서는 안됨 (HARKing, Hypothesis After the Results are Known)

• 탐색형 연구: 전체를 탐색적으로 해석하는 연구
• 가설검증형 연구: 가설을 세우고 이를 검증하는 연구

가설검정을 이해할 때 확인할 항목 (전부 잘못된 명제임)

• 유의미한 결과(p<0.05)란 귀무가설이 틀리다는 뜻이다.
• 유의미가 아닌 결과(p>=0.05)란 귀무가설이 옳다, 또는 귀무가설을 채택해야 한다는 의미다.
• p값이 클수록 귀무가설을 지지한다는 증거이다.
• p>=0.05는 효과가 없음을 관찰했다, 또는 효과가 없다는 것을 입증했다는 것을 뜻한다.
• 통계적으로 유의미하다는 것은 과학적으로 매우 중요한 관계를 밝혔다는 것을 뜻한다.
• 통계적으로 유의미하지 않다면 효과크기가 작다는 뜻이다.