R을 이용한 통계 분석 - 특이점 진단

회귀 분석에서 중요한 가정 중에 하나는 자료에 포함된 관측값에 대해 사용된 모형이 적절해야 한다는 것입니다. 그러나 실제 자료에서는 1~2개의 자료들이 모형을 따르지 않는 경우가 있는데요. 이러한 관측값을 특이점 또는 이상점이라고 부릅니다.

다음은 R을 통해 관측값을 분석해 보도록 하겠습니다.

포브스의 1857년 논문에서 발췌한 것으로 물 끓는 온도는 화씨를, 대기 압력의 각각 단위로 기록된 측정치입니다. 포브스는 물 끓는 온도와 대기 압력에 로그를 취한 값 사이에는 선형관계가 존재한다고 하였습니다. 대기 압력에 로그를 취한 값에 단위를 맞추기 위해 100을 곱한 값을 반응 변수로 하여 분석해 보겠습니다.

forbes.txt

데이터를 일단 불러 옵니다.

> forbes <- read.delim("~/Documents/forbes.txt")

> attach(forbes)

# plot으로 온도와 대기압력에 대한 산점도를 그립니다. 

> plot(temp, log100pres)

#이상치들의 값을 확인하기 위해 identity를 이용해서 클릭해 봅니다.

> identify(temp, log100pres)

위와 같은 산점도를 얻을 수 있습니다. 위의 값에서 보면 12번째 관측값은 특이점이 될 수 있다는 것을 시각적으로 알 수 있습니다.

이제 회귀 분석 적합을 하고 잔차 분석을 통해 위의 데이터의 특징을 확인해 보도록 하겠습니다.

#회귀 모형 적합을 합니다.

> forbes.reg <- lm(log100pres~temp, data=forbes)

> residual <- forbes.reg$resid

#내면 스튜던트 잔차

> istudent.resid <- rstandard(forbes.reg)

#외면 스튜던츠 잔차

> estudent.resid <- rstudent(forbes.reg)

#hat 행렬을 구하기

> hat <- hatvalues(forbes.reg)

#컬럼 바인드로 정렬해주기

> resid.analsys <- cbind(residual, istudent.resid, estudent.resid, hat)

#소수점 4번째자리까지 표시

> resid.analsys <- round(resid.analsys, 4)

> resid.analsys

> sink()

결과는 아래와 같이 나온다.

보면, 12번째 관측값은 잔차가 크기 때문에 특이점이 될 수 있다는 사실을 알았다. 

다음은 car(Companion to Applied Regression) 패키지를 통해서 파악한 이상치의 그래프이다. 

> library(car)

> qq.plot(resid.analsys)

위 자료를 통해서도 12번째 값이 특이점이라는 사실을 알 수 있다.

다음은

outlierTest 패키지를 통해 계산해본 결과이다.

> outlierTest(forbes.reg)

 rstudent unadjusted p-value Bonferonni p

12 12.40691         6.0892e-09   1.0352e-07

위의 결과를 해석하자면 H0 : 12번째 데이터는 outlier가 아니다 가 된다. 하지만, Bonferonni P값이 1.0352e-07이 유의수준 (0.05)를 넘지 않기 때문에 H0을 기각한다. 때문에 12번째 데이터는 특이값으로 판정할 수 있게 된다.

특이점에 대해서 식별할 수 있는 방법은 다양하다. 문제는 이제 이 특이점을 찾은 후에 어떻게 할 것인가에 대해서는 문제에 따라 달라진다. 어떤 경우에는 특이점이 연구의 대상이 될 수없다고 해서 버려질 수도 있고, 어떤 경우에서는 특이점이 연구의 대상일 수도 있다. 특이 점에 대한 원인의 연구를 규명하는 것이 또 다른 문제해결 능력을 요할 수도 있다. 원인을 제거 시켜도 특이점이 없어지지 않는데 임의로 제거 한다는 것은 바람직하지 못하다고 할 수 있다.