Python_Statistics_Regression

선형회귀

• 머신러닝의 회귀는 예측의 성공 확률을 높이는데 목적
• 통계방법의 회귀는 정해진 분포나 가정을 통해 실패 확률을 줄이고 원인을 찾는데 목적

[개념]

• 하나 혹은 그 이상의 원인이 종속변수에 미치는 영향도를 추적
• 변수들 사이의 상관관계를 밝힘
• 관심있는 변수를 예측하거나 추론하기 위해 사용하는 분석 방법

[단순회귀분석의 평가]

• $ SSE $: $ (관측값 - 예측값)^2 $ > 설명되지 않는 변동(잔차)
• $ SSR $: $ (예측값 - 기대값) ^2 $ > 설명할 수 있는 변동
• $ SST $: 총 변동

[출처] https://bioinformaticsandme.tistory.com/70

 

• F-통계량은 설명할 수 있는 변동이 설명하지 못하는 변동보다 클때 값이 크게 나타난다.
• 이는 회귀선이 데이터에 잘 적합했음을 의미한다.
• 따라서 F-통계량이 클수록 귀무가설을 기각할 확률이 올라간다.

 

• 모든 데이터위에 회귀선이 존재하면 $ R^2 $ = 1
• 추정된 회귀선의 기울기가 0이면 $ R^2 $ = 0
• $ R^2 $은 자유도에 영향을 받기 때문에 절대적인 값은 아님
• $ R^2 = r^2 $

[다중회귀분석의 평가]

• 다중회귀분석은 독립변수가 2개 이상일때 사용
• 독립변수가 늘어나면 $ R^2 $ 값이 늘어나게 됨
• 따라서 모델이 독립변수의 수에 따라 증가하게 되는 것에 대한 패널티를 적용
• $ AdjustedR^2  = 1-  \frac {SSR/(n-k-1)}{SST/(n-1)} $

[회귀분석의 검토사항]

1. 모형이 데이터를 잘 적합하고 있는가? > 모형의 잔차가 특정 패턴을 이루고 있지 않아야 한다
   • 잔차가 랜덤하게 분포되어 있는 경우: 등분산성을 만족함
   • 등분산성이 벗어난경우: 분산이 증가하면 가중회귀 직선이 적합
   • 선형성을 벗어나는 경우: 이차항을 추가한 다항 회귀가 적합
   • 잔차가 선형을 이루는 경우: 새로운 독립변수를 추가하는 것이 적합
2. 회귀모형이 통계적으로 유의한가?
3. 모형은 데이터를 얼마나 설명할 수 있는가?
4. 모형의 회귀계수는 유의한가?

[단순선형회귀분석 예시]

• kc_house_data 활용
• 독립변수: sqft_living(평방피트) / 종속변수: price(집값)
• 독립변수가 종속변수에 미치는 영향도 확인

statsmodels.formula.api.ols(parameters)

• 데이터정의: house변수로 할당(21613, 18)
• 필요데이터 추출

#필요데이터 추출
house = house[["price", "sqft_living"]]
print(house.shape)
house.head()

• 독립변수와 종속변수의 선형성 확인

#독립변수와 종속변수간 선형성 확인: 0.702
house.corr()

 

• 모형 적합 및 예측선 시각화

#독립변수와 종속변수 정의
y = house["price"]
x = house[["sqft_living"]]

#모형적합
lr = ols("price ~ sqft_living", data=house).fit()
y_pred = lr.predict(x)

plt.scatter(x, y)
plt.plot(x, y_pred, c="red")
plt.xlabel("sqft_livig")
plt.ylabel("price", rotation=90)
plt.title("Linear Regression Result")
plt.show();

• 회귀분석 검토
  • 모형이 데이터를 잘 적합하고 있는가?
    • 직관적으로 모형(회귀선)이 데이터를 잘 설명하지 못하는 것으로 확인된다.
    • 오차또한 점점 분산이 커지는 특정한 패턴을 이루고 있다.
  • 회귀모형이 적합한가?
    • 귀무가설: 회귀모형이 적합하지 않다.
    • 대립가설: 회귀모형이 적합하다.
    • p-value가 0으로 유의수준 0.5 이하에서 귀무가설을 기각한다.
    • 따라서 모델은 유의하다.
  • 모형은 데이터를 얼마나 설명할 수 있는가?
    • $ R^2 = 0.493 $ 으로 모델은 전체 데이터의 49.3%를 설명할 수 있다.
    • 이는 모델의 성능이 좋지 않음을 나타낸다.
    • 하지만 산업에 따라 해당 수치는 좋을 수도, 좋지 않을수도 있다.
  • 모형 내의 회귀계수는 유의한가?
    • sqft_living변수는 p-value가 0 이므로 .05 이하에서 유의하다.
  • 따라서 $ Price = sqft\_living * 280.6236 - 43580.743094 $ 라는 식이 도출된다.

 

[다중회귀분석]

• 다중회귀분석은 독립변수가 2개 이상
• 다중회귀분석에서는 필수적으로 다중공선성을 확인하고 이를 해결해야 함
• 다중공선성: 독립변수들간의 강한 상관관계가 나타나는 것을 의미
• 다중공선성이 있으면 정확한 회귀계수의 추정이 어려움
• 회귀분석에서는 독립변수의 수가 증가할 수록 모델의 정확도가 올라가는 문제가 발생
• 따라서 다중공선성이 존재한다면 하나의 변수를 제거해 주거나 해당 변수에 패널티를 주어 모델에 영향을 미치는 영향도를 줄여야 함
• 다중공선성을 파악하는 방법
  • 독립변수들간의 상관계수를 구하여 상관성이 0.9 이상인 변수 탐색
  • 다중공선성이 의심되는 두 변수의 회귀분석으로 허용 오차($ 1-R^2 $)를 구했을때, 0.1이하이면 다중공선성 문제가 심각하다고 할 수 있음
  • VIF값이 10 이상이라면 다중공선성이 존재할 것으로 예상 $ VIF = \frac {1}{(1-R^2)} $
  • 변수선택법
    • 종속변수에 영향을 미치는 유의미한 독립변수만을 선택하여 최적의 회귀 방정식을 도출하는 과정이 필요
    • 변수를 선택할때는 모델의 유의성 판단의 근거로 삼았던 F통계량이나 AIC를 기준으로 선택한다.
    • $ AIC = -2ln(L) + 2k $
      • 2ln(L): 모형의 적합도(L: Likelihood function을 의미하고 AIC값이 낮다는 것은 모형의 적합도가 높다는 것을 의미)
      • 2k: 추정된 파라미터의 개수로 해당 모형에 패널티를 주기 위해 사용
    • 단계적 선택법이 정확도가 가장 높음

변수선택법 구분	내용
전진선택법	변수를 하나씩 추가해 가며 모델의 정확도를 높이는 방법
후진제거법	모든 변수를 추가한 상태에서 유의하지 않은 변수들을 제거해 가며 모델의 정확도를 높이는 방법
단계적 선택법	변수 추가, 제거를 반복하며 모델의 성능을 높이는 방법

 

예시) Cars93 데이터 활용 다중회귀분석 실시

• 데이터정의: Cars93(93, 27)
• 독립변수와 종속변수 정의
  • 독립변수: EngineSize, RPM, Weight, Length, MPG.city, MPG.highway
  • 종속변수: Price

#필요데이터 추출
cars = cars[["EngineSize", "RPM", "MPG.city", "MPG.highway", "Price", "Weight", "Length"]]
print(cars.shape)
cars.head()

#컬럼명 전처리
cars.columns = cars.columns.str.replace(".", "")
cars.columns = cars.columns.str.lower()
cars.head(3)

#독립변수와 종속변수 분리
y = cars["price"]
x_cols = cars.columns.difference(["price"])
x = cars[x_cols]

 

• 모델 적합
  • 기본적인 다중회귀분석을 실시한 결과
  • 모델의 p-value .05 하에서 유의
  • EngineSize, RPM 두개의 변수만 .05 하에서 유의
  • $ Adj.R^2 $ 가 54.2%로 나타남
  • 다중공선성을 확인하고 처리한 결과를 확인

#모델적합
from statsmodels.formula.api import ols

lr = ols(f"price ~ {'+'.join(x_cols)}", data=cars).fit()
print(lr.summary())

 

• 다중공선성 확인을 위해 변수간 상관계수 확인
  • Mpgcity와 Mpghighway 두 변수의 상관관계가 0.9로 매우 높게 나타남
  • 둘중 한 변수 제거해야 함

 

• vif를 통해 어떤 변수를 제거할 것인지 확인
  • from patsy import dmatrices
  • from statsmodels.stats.outliers_influence import variance_inflation_factor
  • MPGcity의 변수가 가장 높게 나타남으로 제거해야 함을 알 수 있음

#vif 산출
from patsy import dmatrices
from statsmodels.stats.outliers_influence import variance_inflation_factor

y, X = dmatrices(f"price ~ {'+'.join(x_cols)}", data=cars, return_type="dataframe")

X

vif_list = []
for i in range(1, len(X.columns)):
  vif_list.append(
      [variance_inflation_factor(X.values, i), X.columns[i]]
  )

vif_df = pd.DataFrame(vif_list, columns=["vif", "var"])
vif_df

dmatrices > X / vif산출

 

• 변수를 제거하고 다시 다중선형회귀분석 실시
  • $ Adj.R^2 $ 값이 눈에띄게 나아지지는 않음
  • $ AIC $ 값도 눈에 띄게 나아지지 않음
  • 상관계수가 높았던 MPGhighway 의 p-value값이 낮아지고 나머지 변수들의 p-value값도 변화하는 것을 볼 수 있음
  • 이와같이 다중공선성을 해결하지 않으면 유의한 변수임에도 불구하고 유의하지 않다고 여겨질 가능성, 유의하지 않은 변수임에도 유의하다고 여겨질 가능성이 다분함으로 꼭 해결해야 함

#다시 다중 선형회귀분석 진행
from statsmodels.formula.api import ols

lr = ols(f"price ~ enginesize + length + mpghighway + rpm + weight", data=cars).fit()
print(lr.summary())

좌: 변수제거 전 / 우: 변수제거 후

 

[변수선택법 실시]

• 전진선택법: 상수에서 시작해 변수를 하나씩 늘려가며 AIC를 계산하는 방법
  • 현 데이터에서 전진선택법의 결과도 모든 변수를 활용해 분석한 회귀분석의 AIC와 동일하게 나왔음을 알 수 있다.

#필요데이터 추출
cars = cars[["enginesize", "rpm", "mpgcity", "mpghighway", "price", "weight", "length"]]
print(cars.shape)
cars.head()

#독립변수와 종속변수 분리
y = cars["price"]
x_cols = cars.columns.difference(["price"])
x = cars[x_cols]

print(y.shape, x.shape)

#상수항 붙도록 추가
from patsy import dmatrices

y, X = dmatrices(f"price ~ {'+'.join(cars.columns.difference(['price']))}", data=cars, return_type="dataframe")
print(X.shape)

X.head()

 

#전체변수로 진행
import statsmodels.api as sm

lr = sm.OLS(y, X).fit()
aic = lr.aic
result = {"model": lr.model.exog_names, "AIC": aic}
result_df = pd.DataFrame([result])
result_df

#전진선택법 함수 정의
def forward(y, X, vars):
  features = [f for f in X.columns.difference(["Intercept"]) if f not in vars]
  result_list = []
  for i in range(len(features)):
    lr = sm.OLS(y, X[vars + [features[i]] + ["Intercept"]]).fit()
    aic = lr.aic
    result = {"model": lr.model.exog_names, "AIC": aic}
    result_list.append(result)

  result_df = pd.DataFrame([result])

  best_model = result_df.loc[result_df["AIC"].argmin()]

  return best_model

#----------전진선택법: 상수에서 시작해서 모든 변수 하나로 합치기
#-----상수만
lr = sm.OLS(y, X["Intercept"]).fit()
aic = lr.aic
results = {"model": lr.model.exog_names, "AIC": aic}
print(results)

#----전진선택
forward_df = pd.DataFrame(columns=["model", "AIC"])
vars = []
for i in range(len(X.columns.difference(["Intercept"]))):
  print()
  best_model = forward(y, X, vars)
  print(best_model)

  forward_df.loc[i] = best_model
  vars = best_model.model #weight, Intercept
  
  vars = [v for v in vars if v != "Intercept"] #weight
  print(vars)
forward_df

 

• 후진제거법: 모든 변수에서 하나씩 줄여가며 AIC를 확인하는 방법
  • enginesize, rmp, weight, Intercept로 분석한 AIC의 결과가 619.09765로 가장 낮게 나타남을 알 수 있다.

def backward(X, y, vars):
  result_list = []
  for combo in itertools.combinations(vars, len(vars)-1):
    lr = sm.OLS(y, X[list(combo) + ["Intercept"]]).fit()
    aic = lr.aic
    result = {"model": lr.model.exog_names, "AIC": aic}
    result_list.append(result)
  result_df = pd.DataFrame(result_list)
  best_model = result_df.loc[result_df["AIC"].argmin()]

  return best_model

#후진제거법
#일단 모든 변수로 돌려서 aic 점수 확인한다.
full_lr = sm.OLS(y, X).fit() #상수와 6개 변수 = 7
full_aic = full_lr.aic
full_result = {"model": full_lr.model.exog_names, "AIC": full_aic}
print(full_result)

vars = X.columns.difference(["Intercept"]) # 6개

#결과 저장할 데이터 프레임 생성
best_df = pd.DataFrame(
    columns=["model", "AIC"],
    index=range(1, len(vars)+1)
)
best_df.loc[len(vars)] = full_result #모든 변수로 분석한 결과 저장

while len(vars) > 1:
  best_result = backward(X, y, vars) #후진제거법 실시
  print()
  print(best_result)

  if best_result["AIC"] > full_aic:
    break

  best_df.loc[len(vars)-1] = best_result
  full_aic = best_result["AIC"]

  vars = best_df.loc[len(vars)-1, "model"]
  vars = [v for v in vars if v != "Intercept"]

best_df

 

• 단계선택법: 전진과 후진을 AIC기준으로 비교해가며 실행

# 단계적 선택법
# 전진선택법 부터 시작

vars = [] 
result_before = sm.OLS(y, X[vars + ["Intercept"]]).fit() #intercept만 가지고 실시
print(f"result_before: {result_before.model.exog_names}, {result_before.aic}")


for i in range(len(X.columns.difference(["Intercept"]))): #6번 돌기
  print()
  print()
  print("[",i, c,"]")
  forward_result = forward(y, X, vars) #전진선택법
  print(f"*****forward_Result*****\n{forward_result}")
  print("forward!!!!")

  vars = forward_result["model"] #vars업데이트
  vars = [v for v in vars if v != "Intercept"]
  print()
  print(f"updated_vars: {vars}")
  print()

  backward_result = backward(X, y, vars) #후진제거법
  print(f"*****backward_result*****\n{backward_result}")

  if backward_result["AIC"] < forward_result["AIC"]:    
    result_before.aic = backward_result.AIC
    print(f"result_before: {result_before.aic}")
    vars = backward_result["model"]
    vars = [v for v in vars if v != "Intercept"]
    print("backward!!!!!")

  
  if backward_result["AIC"] > result_before.aic:
    break
    print("backward_result is better than result_before")
  else:
    result_before.aic = backward_result.AIC
    print(f"updated result_before: {result_before.aic}")