Я заинтересован в определении количества значимых паттернов, вытекающих из анализа основных компонентов (PCA) или анализа эмпирических ортогональных функций (EOF). Я особенно заинтересован в применении этого метода к климатическим данным. Поле данных представляет собой матрицу MxN, где М - это измерение времени (например, дни), а N - пространственное измерение (например, положения в долготе / широте). Я читал о возможном методе начальной загрузки для определения значимых ПК, но не смог найти более подробное описание. До сих пор я применял «Правило большого пальца Норта» (North et al ., 1982), чтобы определить это ограничение, но мне было интересно, доступен ли более надежный метод.
В качестве примера:
###Generate data
x <- -10:10
y <- -10:10
grd <- expand.grid(x=x, y=y)
#3 spatial patterns
sp1 <- grd$x^3+grd$y^2
tmp1 <- matrix(sp1, length(x), length(y))
image(x,y,tmp1)
sp2 <- grd$x^2+grd$y^2
tmp2 <- matrix(sp2, length(x), length(y))
image(x,y,tmp2)
sp3 <- 10*grd$y
tmp3 <- matrix(sp3, length(x), length(y))
image(x,y,tmp3)
#3 respective temporal patterns
T <- 1:1000
tp1 <- scale(sin(seq(0,5*pi,,length(T))))
plot(tp1, t="l")
tp2 <- scale(sin(seq(0,3*pi,,length(T))) + cos(seq(1,6*pi,,length(T))))
plot(tp2, t="l")
tp3 <- scale(sin(seq(0,pi,,length(T))) - 0.2*cos(seq(1,10*pi,,length(T))))
plot(tp3, t="l")
#make data field - time series for each spatial grid (spatial pattern multiplied by temporal pattern plus error)
set.seed(1)
F <- as.matrix(tp1) %*% t(as.matrix(sp1)) +
as.matrix(tp2) %*% t(as.matrix(sp2)) +
as.matrix(tp3) %*% t(as.matrix(sp3)) +
matrix(rnorm(length(T)*dim(grd)[1], mean=0, sd=200), nrow=length(T), ncol=dim(grd)[1]) # error term
dim(F)
image(F)
###Empirical Orthogonal Function (EOF) Analysis
#scale field
Fsc <- scale(F, center=TRUE, scale=FALSE)
#make covariance matrix
C <- cov(Fsc)
image(C)
#Eigen decomposition
E <- eigen(C)
#EOFs (U) and associated Lambda (L)
U <- E$vectors
L <- E$values
#projection of data onto EOFs (U) to derive principle components (A)
A <- Fsc %*% U
dim(U)
dim(A)
#plot of top 10 Lambda
plot(L[1:10], log="y")
#plot of explained variance (explvar, %) by each EOF
explvar <- L/sum(L) * 100
plot(explvar[1:20], log="y")
#plot original patterns versus those identified by EOF
layout(matrix(1:12, nrow=4, ncol=3, byrow=TRUE), widths=c(1,1,1), heights=c(1,0.5,1,0.5))
layout.show(12)
par(mar=c(4,4,3,1))
image(tmp1, main="pattern 1")
image(tmp2, main="pattern 2")
image(tmp3, main="pattern 3")
par(mar=c(4,4,0,1))
plot(T, tp1, t="l", xlab="", ylab="")
plot(T, tp2, t="l", xlab="", ylab="")
plot(T, tp3, t="l", xlab="", ylab="")
par(mar=c(4,4,3,1))
image(matrix(U[,1], length(x), length(y)), main="eof 1")
image(matrix(U[,2], length(x), length(y)), main="eof 2")
image(matrix(U[,3], length(x), length(y)), main="eof 3")
par(mar=c(4,4,0,1))
plot(T, A[,1], t="l", xlab="", ylab="")
plot(T, A[,2], t="l", xlab="", ylab="")
plot(T, A[,3], t="l", xlab="", ylab="")
И вот метод, который я использовал для определения значимости ПК. Основное правило заключается в том, что разница между соседними лямбдами должна быть больше, чем связанная с ними ошибка.
###Determine significant EOFs
#North's Rule of Thumb
Lambda_err <- sqrt(2/dim(F)[2])*L
upper.lim <- L+Lambda_err
lower.lim <- L-Lambda_err
NORTHok=0*L
for(i in seq(L)){
Lambdas <- L
Lambdas[i] <- NaN
nearest <- which.min(abs(L[i]-Lambdas))
if(nearest > i){
if(lower.lim[i] > upper.lim[nearest]) NORTHok[i] <- 1
}
if(nearest < i){
if(upper.lim[i] < lower.lim[nearest]) NORTHok[i] <- 1
}
}
n_sig <- min(which(NORTHok==0))-1
plot(L[1:10],log="y", ylab="Lambda (dots) and error (vertical lines)", xlab="EOF")
segments(x0=seq(L), y0=L-Lambda_err, x1=seq(L), y1=L+Lambda_err)
abline(v=n_sig+0.5, col=2, lty=2)
text(x=n_sig, y=mean(L[1:10]), labels="North's Rule of Thumb", srt=90, col=2)
Я посчитал полезным раздел главы Björnsson and Venegas ( 1997 ) о значимых тестах - они относятся к трем категориям тестов, из которых, вероятно, я надеюсь использовать доминирующий тип дисперсии . Относится к типу подхода Монте-Карло, в котором тасуется измерение времени и пересчитывается лямбда-многовариант по многим перестановкам. von Storch и Zweiers (1999) также ссылаются на тест, который сравнивает лямбда-спектр с эталонным спектром «шума». В обоих случаях я немного не уверен в том, как это можно сделать, а также в том, как проводится тест значимости с учетом доверительных интервалов, определенных перестановками.
Спасибо за вашу помощь.
Ссылки: Бьорнссон, Х. и Венегас, С.А. (1997). «Руководство по EOF и SVD анализу климатических данных», Университет Макгилла, Отчет CCGCR № 97-1, Montréal, Québec, 52pp. http://andvari.vedur.is/%7Efolk/halldor/PICKUP/eof.pdf
GR North, TL Bell, RF Cahalan и FJ Moeng. (1982). Ошибки выборки при оценке эмпирических ортогональных функций. ПН Wea. Rev., 110: 699–706.
фон Шторх, H, Цвайерс, FW (1999). Статистический анализ в исследовании климата. Издательство Кембриджского университета.
источник
Ответы:
Я собираюсь немного продвинуть диалог, хотя это мой вопрос. Прошло 6 месяцев с тех пор, как я спросил об этом, и, к сожалению, не было дано полных ответов, я постараюсь обобщить то, что я собрал до сих пор, и посмотрим, сможет ли кто-нибудь уточнить оставшиеся вопросы. Пожалуйста, извините за длинный ответ, но я не вижу другого пути ...
Сначала я продемонстрирую несколько подходов, использующих возможно лучший синтетический набор данных. Это взято из статьи Беккерса и Риксона ( 2003 ), иллюстрирующей использование алгоритма для проведения EOF на gappy данных. Я воспроизвел алгоритм на R, если кому-то интересно ( ссылка ).
Синтетический набор данных:
Итак, поле истинных данных
Xt
состоит из 9 сигналов, и я добавил к нему некоторый шум, чтобы создать наблюдаемое полеXp
, которое будет использовано в примерах ниже. EOF определяются как таковые:EOF
Следуя примеру, который я использовал в своем первоначальном примере, я буду определять «значимые» EOFs по эмпирическому правилу.
Правило Севера
Поскольку лямбда-значения 2: 4 очень близки друг к другу по амплитуде, они считаются незначительными по эмпирическому правилу - то есть их соответствующие EOF-шаблоны могут перекрываться и смешиваться, учитывая их аналогичные амплитуды. Это прискорбно, учитывая, что мы знаем, что в поле действительно существует 9 сигналов.
Более субъективный подход состоит в том, чтобы просмотреть преобразованные в лог значения Lambda («график Scree») и затем подогнать регрессию к конечным значениям. Затем можно визуально определить, на каком уровне лямбда-значения лежат выше этой линии:
Сюжетный сюжет
Итак, 5 ведущих EOF лежат выше этой линии. Я пробовал этот пример, когда
Xp
не было добавлено никакого дополнительного шума, и результаты показывают все 9 исходных сигналов. Итак, незначительность EOFs 6: 9 связана с тем, что их амплитуда ниже, чем шум в поле.Более объективный метод - критерии «Правила N» по Overland и Preisendorfer (1982). В
wq
пакете есть реализация , которую я покажу ниже.Правило N
Правило N идентифицировало 4 значимых EOF. Лично мне нужно лучше понять этот метод; Почему можно измерить уровень ошибки на основе случайного поля, которое не использует такое же распределение, как в
Xp
? Одним из вариантов этого метода будет повторная выборка данных,Xp
чтобы каждый столбец переставлялся случайным образом. Таким образом, мы гарантируем, что полная дисперсия случайного поля равнаXp
. Повторная выборка много раз позволяет нам вычислить базовую ошибку разложения.Монте-Карло со случайным полем (т.е. сравнение нулевой модели)
Опять же, 4 EOF выше распределений для случайных полей. Я обеспокоен этим подходом и правилом N в том, что они не относятся к доверительным интервалам значений лямбда; Например, высокое первое лямбда-значение автоматически приведет к меньшему количеству отклонений, которые будут объясняться концевыми. Таким образом, лямбда, рассчитанная по случайным полям, всегда будет иметь менее крутой наклон и может привести к выбору слишком небольшого количества значимых EOF. [ПРИМЕЧАНИЕ:
eofNum()
функция предполагает, что EOF рассчитываются по корреляционной матрице. Это число может отличаться, если используется, например, ковариационная матрица (центрированные, но не масштабированные данные).]Наконец, @ tomasz74 упомянул статью Babamoradi et al. (2013), который я кратко рассмотрел. Это очень интересно, но, кажется, больше сосредоточено на расчете КИ нагрузок и коэффициентов EOF, а не на лямбду. Тем не менее, я считаю, что это может быть принято для оценки лямбда-ошибки с использованием той же методологии. Начальная выборка выполняется для поля данных путем повторной выборки строк, пока не будет создано новое поле. Одну и ту же строку можно пересчитать более одного раза, что является непараметрическим подходом, и нет необходимости делать предположения о распределении данных.
Начальная загрузка лямбда-значений
Хотя это может быть более надежным, чем эмпирическое правило Норда, для расчета погрешности значений лямбда, теперь я считаю, что вопрос о значении EOF сводится к различным мнениям о том, что это значит. Для эмпирических правил и методов начальной загрузки Севера значение, по-видимому, больше зависит от того, перекрывается ли значение лямбда или нет. Если есть, то эти EOF могут смешиваться в своих сигналах и не представлять «истинные» шаблоны. С другой стороны, эти два EOF могут описывать значительную дисперсию (по сравнению с разложением случайного поля - например, Правило N). Поэтому, если кто-то заинтересован в фильтрации шума (т. Е. Через усечение EOF), тогда будет достаточно правило N. Если кто-то заинтересован в определении реальных шаблонов в наборе данных, то более строгие критерии лямбда-перекрытия могут быть более надежными.
Опять же, я не эксперт в этих вопросах, поэтому я все еще надеюсь, что кто-то с большим опытом может добавить к этому объяснению.
Ссылки:
Беккерс, Жан-Мари и М. Риксен. «Расчеты EOF и заполнение данных из неполных наборов океанографических данных». Журнал атмосферных и океанических технологий 20.12 (2003): 1839-1856.
Overland, J., R. Preisendorfer, Тест значимости для основных компонентов, применяемых в циклонической климатологии, пн. Wea. Rev., 110, 1-4, 1982.
источник