기본 콘텐츠로 건너뛰기

[Linear Algebra] 직교적 투영(Orthogonal Projection)

[data analysis] Hat 행렬

Hat 행렬

관련된 내용

hat 행렬(H)은 회귀모델을 적용하여 $\hat{y}$을 추정하는 과정에서 설명변수의 영향을 나타내기 위해 작성합니다. 식 1과 같이 최소제곱 방법에 의한 회귀계수를 적용한 추정치의 추정 과정에서 설명변수만으로 구성된 부분을 hat 행렬(P)이라 합니다.

\begin{align}\hat{y} & = Xb\\&=X(X^TX)^{-1}X^Ty\\& = Py\\ H&= X(X^TX)^{-1}X^T\\& = P\end{align} (식 1)

이상치는 그 데이터가 보이는 일반적인 경향을 따르지 않는 값(들)을 의미합니다. 회귀분석시 일반적으로 이상치는 반응변수 값의 측면에서 고려되지만 설명변수 측면에서 극단의 값(x값)을 가진다면 그 데이터 지점은 높은 레버리지(leverage)를 가진다고 말합니다. 이 레버리지 즉, 각 샘플의 설명변수의 영향은 hat 행렬의 대각요소의 값으로 나타낼 수 있습니다.

hat 행렬은 이상치를 판단하기 위한 기본 값들로 이를 기준으로 여러 지표를 계산할 수 있습니다. statsmodel.api.OLS()에 의해 생성되는 모델의 get_influence() 메소드는 이상치 판단을 위한 여러 지표들의 결과를 나타냅니다. 이 결과는 .summary_frame() 메소드로 확인할 수 있습니다.

기사 회귀모형에서 이상치(outlier) 파악편에서 구축한 회귀모델 reg에 대해 이 메소드들을 실행하면 다음과 같습니다.

import numpy as np 
import pandas as pd 
from sklearn.preprocessing import StandardScaler 
import matplotlib.pyplot as plt
import FinanceDataReader as fdr
import statsmodels.api as sm
from scipy import stats
st=pd.Timestamp(2021,1, 1)
et=pd.Timestamp(2024, 5, 10)
kos=fdr.DataReader('KS11',st, et)[["Open","Close"]]
kos.index=range(len(kos))
X=kos.values[:,0].reshape(-1,1)
y=kos.values[:,1].reshape(-1,1)

#독립변수 정규화(표준화)
xScaler=StandardScaler().fit(X)
X_n=xScaler.transform(X)
#반응변수 정규화(표준화)
yScaler=StandardScaler().fit(y)
y_n=yScaler.transform(y)

X_n0=sm.add_constant(X_n)
X_n0.shape, y_n.shape
reg=sm.OLS(y_n, X_n0).fit()
influence=reg.get_influence()
infSummary=influence.summary_frame()
infSummary.round(4).head(3)
dfb_const dfb_x1 cooks_d standard_resid hat_diag dffits_internal student_resid dffits
0 0.1106 0.0620 0.0079 3.1599 0.0016 0.1261 3.1773 0.1268
1 0.0755 0.0597 0.0046 2.1649 0.0020 0.0961 2.1698 0.0963
2 -0.0337 -0.0323 0.0011 -0.9695 0.0023 -0.0467 -0.9694 -0.0467

이 블로그의 인기 게시물

유사변환과 대각화

내용 유사변환 유사행렬의 특성 대각화(Diagonalization) 유사변환(Similarity transformation) 유사변환 n×n 차원의 정방 행렬 A, B 그리고 가역 행렬 P 사이에 식 1의 관계가 성립하면 행렬 A와 B는 유사하다고 하며 이 변환을 유사 변환 (similarity transformation)이라고 합니다. $$\begin{equation}\tag{1} A = PBP^{-1} \Leftrightarrow P^{-1}AP = B \end{equation}$$ 식 1의 유사 변환은 다음과 같이 고유값을 적용하여 특성 방정식 형태로 정리할 수 있습니다. $$\begin{align} B - \lambda I &= P^{-1}AP – \lambda P^{-1}P\\ &= P^{-1}(AP – \lambda P)\\ &= P^{-1}(A - \lambda I)P \end{align}$$ 위 식의 행렬식은 다음과 같이 정리됩니다. $$\begin{align} &\begin{aligned}\textsf{det}(B - \lambda I ) & = \textsf{det}(P^{-1}(AP – \lambda P))\\ &= \textsf{det}(P^{-1}) \textsf{det}((A – \lambda I)) \textsf{det}(P)\\ &= \textsf{det}(P^{-1}) \textsf{det}(P) \textsf{det}((A – \lambda I))\\ &= \textsf{det}(A – \lambda I)\end{aligned}\\ &\begin{aligned}\because \; \textsf{det}(P^{-1}) \textsf{det}(P) &= \textsf{det}(P^{-1}P)\\ &= \t

[matplotlib] 히스토그램(Histogram)

히스토그램(Histogram) 히스토그램은 확률분포의 그래픽적인 표현이며 막대그래프의 종류입니다. 이 그래프가 확률분포와 관계가 있으므로 통계적 요소를 나타내기 위해 많이 사용됩니다. plt.hist(X, bins=10)함수를 사용합니다. x=np.random.randn(1000) plt.hist(x, 10) plt.show() 위 그래프의 y축은 각 구간에 해당하는 갯수이다. 빈도수 대신 확률밀도를 나타내기 위해서는 위 함수의 매개변수 normed=True로 조정하여 나타낼 수 있다. 또한 매개변수 bins의 인수를 숫자로 전달할 수 있지만 리스트 객체로 지정할 수 있다. 막대그래프의 경우와 마찬가지로 각 막대의 폭은 매개변수 width에 의해 조정된다. y=np.linspace(min(x)-1, max(x)+1, 10) y array([-4.48810153, -3.54351935, -2.59893717, -1.65435499, -0.70977282, 0.23480936, 1.17939154, 2.12397372, 3.0685559 , 4.01313807]) plt.hist(x, y, normed=True) plt.show()

R 미분과 적분

내용 expression 미분 2차 미분 mosaic를 사용한 미분 적분 미분과 적분 R에서의 미분과 적분 함수는 expression()함수에 의해 생성된 표현식을 대상으로 합니다. expression expression(문자, 또는 식) 이 표현식의 평가는 eval() 함수에 의해 실행됩니다. > ex1<-expression(1+0:9) > ex1 expression(1 + 0:9) > eval(ex1) [1] 1 2 3 4 5 6 7 8 9 10 > ex2<-expression(u, 2, u+0:9) > ex2 expression(u, 2, u + 0:9) > ex2[1] expression(u) > ex2[2] expression(2) > ex2[3] expression(u + 0:9) > u<-0.9 > eval(ex2[3]) [1] 0.9 1.9 2.9 3.9 4.9 5.9 6.9 7.9 8.9 9.9 미분 D(표현식, 미분 변수) 함수로 미분을 실행합니다. 이 함수의 표현식은 expression() 함수로 생성된 객체이며 미분 변수는 다음 식의 분모의 변수를 의미합니다. $$\frac{d}{d \text{변수}}\text{표현식}$$ 이 함수는 어떤 함수의 미분의 결과를 표현식으로 반환합니다. > D(expression(2*x^3), "x") 2 * (3 * x^2) > eq<-expression(log(x)) > eq expression(log(x)) > D(eq, "x") 1/x > eq2<-expression(a/(1+b*exp(-d*x))); eq2 expression(a/(1 + b * exp(-d * x))) > D(eq2, "x") a * (b * (exp(-d * x) * d))/(1 + b