2016-01-19 29 views
1

我正在從Michael Faraway的線性模型(其中R(第11章,第160頁))開始PCA部分的工作。使用QQ Plot結果刪除現有的新異常值

PCA分析對異常值敏感,Mahalanobis距離有助於我們識別它們。

作者通過繪製馬哈拉諾比斯距離和卡方分佈的分位數來檢查異常值。

if require(faraway)==F install.packages("faraway"); require(faraway) 
data(fat, package='faraway') 
cfat <- fat[,9:18] 

n <- nrow(cfat); p <- ncol(cfat) 
plot(qchisq(1:n/(n+1),p), sort(md), xlab=expression(paste(chi^2, 
                  "quantiles")), 
ylab = "Sorted Mahalanobis distances") 
abline(0,1) 

我找準穴位:

identify(qchisq(1:n/(n+1),p), sort(md)) 

看來,異常值是行242:252。我刪除這些離羣值,重新創建QQ圖:

cfat.mod <- cfat[-c(242:252),] #remove outliers 
robfat <- cov.rob(cfat.mod) 
md <- mahalanobis(cfat.mod, center=robfat$center, cov=robfat$cov) 
n <- nrow(cfat.mod); p <- ncol(cfat.mod) 
plot(qchisq(1:n/(n+1),p), sort(md), xlab=expression(paste(chi^2, 
                  "quantiles")), 
    ylab = "Sorted Mahalanobis distances") 
abline(0,1) 

identify(qchisq(1:n/(n+1),p), sort(md)) 

唉,現在看來,一個新的點的集合(行234:241),是現在異常。每次我刪除更多的異常值時,都會發生這種情況。

期待理解我做錯了什麼。

+0

我們需要'md'解決。 –

+1

你可能沒有做錯任何事情。這是錯誤指定模型的常見問題。 –

+0

同意@Alex,您將根據排除舊異常值的「新」數據集重新計算異常值。 –

回答

2

要正確識別點,請確保標籤對應於數據中點的位置。函數ordersortindex.return=TRUE將給出排序後的索引。以下是一個例子,任意刪除大於閾值的點數爲md

## Your data 
data(fat, package='faraway') 
cfat <- fat[, 9:18] 
n <- nrow(cfat) 
p <- ncol(cfat) 
md <- sort(mahalanobis(cfat, colMeans(cfat), cov(cfat)), index.return=TRUE) 
xs <- qchisq(1:n/(n+1), p) 
plot(xs, md$x, xlab=expression(paste(chi^2, 'quantiles'))) 

## Use indices in data as labels for interactive identify 
identify(xs, md$x, labels=md$ix) 

## remove those with md>25, for example 
inds <- md$x > 25 
cfat.mod <- cfat[-md$ix[inds], ] 
nn <- nrow(cfat.mod) 
md1 <- mahalanobis(cfat.mod, colMeans(cfat.mod), cov(cfat.mod)) 

## Plot the new data 
par(mfrow=c(1, 2)) 
plot(qchisq(1:nn/(nn+1), p), sort(md1), xlab='chisq quantiles', ylab='') 
abline(0, 1, col='red') 
car::qqPlot(md1, distribution='chisq', df=p, line='robust', main='With car::qqPlot') 

​​