Смещенная начальная загрузка: можно ли центрировать КИ вокруг наблюдаемой статистики?

Это похоже на Bootstrap: оценка находится вне доверительного интервала

У меня есть некоторые данные, которые представляют количество генотипов в популяции. Я хочу оценить генетическое разнообразие, используя индекс Шеннона, а также создать доверительный интервал с помощью начальной загрузки. Я заметил, однако, что оценка с помощью начальной загрузки имеет тенденцию быть чрезвычайно смещенной и приводит к доверительному интервалу, который находится вне моей наблюдаемой статистики.

Ниже приведен пример.

# Shannon's index
H <- function(x){
  x <- x/sum(x)
  x <- -x * log(x, exp(1))
  return(sum(x, na.rm = TRUE))
}
# The version for bootstrapping
H.boot <- function(x, i){
  H(tabulate(x[i]))
}

Генерация данных

set.seed(5000)
X <- rmultinom(1, 100, prob = rep(1, 50))[, 1]

расчет

H(X)

## [1] 3.67948

xi <- rep(1:length(X), X)
H.boot(xi)

## [1] 3.67948

library("boot")
types <- c("norm", "perc", "basic")
(boot.out <- boot::boot(xi, statistic = H.boot, R = 1000L))

## 
## CASE RESAMPLING BOOTSTRAP FOR CENSORED DATA
## 
## 
## Call:
## boot::boot(data = xi, statistic = H.boot, R = 1000)
## 
## 
## Bootstrap Statistics :
##     original     bias    std. error
## t1*  3.67948 -0.2456241  0.06363903

Генерация КИ с коррекцией смещения

boot.ci(boot.out, type = types)

## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = boot.out, type = types)
## 
## Intervals : 
## Level      Normal              Basic              Percentile     
## 95%   ( 3.800,  4.050 )   ( 3.810,  4.051 )   ( 3.308,  3.549 )  
## Calculations and Intervals on Original Scale

Предполагая, что дисперсия t может быть использована для дисперсии t0 .

norm.ci(t0 = boot.out$t0, var.t0 = var(boot.out$t[, 1]))[-1]

## [1] 3.55475 3.80421

Было бы правильно сообщить о КИ с центром в момент времени t0 ? Есть ли лучший способ для создания начальной загрузки?

r confidence-interval bootstrap bias diversity ZNK
источник

Ответы:

В установке, заданной OP, интересующим параметром является энтропия Шеннона которая является функцией вектора вероятности . Оценщик, основанный на выборках ( в моделировании), является оценщиком плагина Образцы были получены с использованием равномерного распределения, для которого энтропия Шеннона равнаТак как энтропия Шеннона максимизируется при равномерном распределении, оценка подключаемого модуля должна быть смещена вниз . Моделирование показывает, что

θ (p) = - \sum_{i = 1}^{50} p_{i} \log p_{i},

$\theta(\mathbf{p}) = - \sum_{i = 1}^{50} p_i \log p_i,$

p \in R^{50}

$\mathbf{p} \in \mathbb{R}^{50}$

n

$n$

n = 100

$n = 100$

{\hat{θ}}_{n} = θ ({\hat{p}}_{n}) = - \sum_{i = 1}^{50} {\hat{p}}_{n, i} \log {\hat{p}}_{n, i} .

$\hat{\theta}_n = \theta(\hat{\mathbf{p}}_n) = - \sum_{i=1}^{50} \hat{p}_{n,i} \log \hat{p}_{n,i}.$

\log (50) = 3.912.

$\log(50) = 3.912.$

b i a s ({\hat{θ}}_{100}) ≃ - 0.28

$\mathrm{bias}(\hat{\theta}_{100}) \simeq -0.28$ тогда как . Оценка подключаемого модуля непротиворечива, но метод не применяется для равномерного распределения , поскольку производная энтропии Шеннона равна 0. Таким образом, для этого конкретного выбора , доверительные интервалы, основанные на асимптотических аргументах, не очевидны.

b i a s ({\hat{θ}}_{500}) ≃ - 0.05

$\mathrm{bias}(\hat{\theta}_{500}) \simeq -0.05$

Δ

$\Delta$

p

$\mathbf{p}$

p

$\mathbf{p}$

Процентный интервал основан на распределении где - это оценка, полученная из выборки наблюдений из . В частности, это интервал от квантиля 2,5% до квантиля 97,5% для распределения . Как показывает имитация начальной загрузки OP, явно также смещена вниз как оценка , что приводит к тому, что процентильный интервал совершенно неправильно. $\theta(\mathbf{p}_n^*)$ $\mathbf{p}_n^*$ $n$ $\hat{\mathbf{p}}_n$ $\theta(\mathbf{p}_n^*)$ $\theta(\mathbf{p}_n^*)$ $\theta(\hat{\mathbf{p}}_n)$

Для основного (и нормального) интервала роли квантилей меняются местами. Это означает, что интервал кажется разумным (он охватывает 3,912), хотя интервалы, выходящие за пределы 3,912, не имеют логического значения. Более того, я не знаю, будет ли базовый интервал иметь правильное покрытие. Его обоснование основано на следующем приблизительном распределительном тождестве:

θ (p_{n}^{*}) - θ ({\hat{p}}_{n}) \overset{D}{≃} θ ({\hat{p}}_{n}) - θ (p),

$\theta(\mathbf{p}_n^*) - \theta(\hat{\mathbf{p}}_n) \overset{\mathcal{D}}{\simeq} \theta(\hat{\mathbf{p}}_n) - \theta(\mathbf{p}),$ что может быть сомнительно для (относительно) малого такого как .

n

$n$

n = 100

$n = 100$

Последнее предложение ОП относительно стандартного интервала на основе ошибок также не будет работать из-за большого смещения. Это может работать для оценки с поправкой на смещение, но в первую очередь вам нужны правильные стандартные ошибки для оценки с поправкой на смещение. $\theta(\hat{\mathbf{p}}_n) \pm 1.96\hat{\mathrm{se}}_n$

Я бы рассмотрел интервал вероятности на основе логарифмической вероятности профиля для . Боюсь, что я не знаю простого способа вычисления вероятности записи профиля для этого примера, за исключением того, что вам нужно максимизировать вероятность записи журнала над для различных фиксированных значений . $\theta(\mathbf{p})$ $\mathbf{p}$ $\theta(\mathbf{p})$

NRH
источник

Проблема смещения с использованием «подключаемого» оценщика для энтропии ценилась десятилетиями. В этой статье анализируются менее предвзятые оценки. Поправка смещения до порядка , которая датируется 1955 годом (см. Уравнение 4 в связанной статье), может быть применена к случаю, представленному ФП. Коррекция составляет 0,245, почти идентична смещению, выявленному бутстрапом. Возможно, здесь следует использовать бутстрап для оценки самой энтропии, а не только ее доверительных границ.

1 / n

$1/n$

EdM

@EdM это очень полезная информация. Я не знал литературы по этой конкретной проблеме смещения. Было бы действительно полезно, если бы вы могли превратить комментарий в ответ, объясняющий исправление смещения и то, как его можно использовать с начальной загрузкой, скажем, для получения доверительных интервалов.

NRH

Я тоже не знал этой литературы, пока не появился этот вопрос и ваш ответ. Что несколько смущает, так как энтропия Шеннона часто используется в качестве меры в моей области биомедицинской науки. Я посмотрю, что я могу собрать в качестве дополнительного ответа.

EdM

Увеличение количества образцов начальной загрузки не поможет. Он должен быть достаточно большой , чтобы можно было надежно оценить величины интереса для распределения , скажем, но в противном случае увеличения числа образцов бутстраповских не удалить смещение или сделать доверие больше подходит.

θ (p_{n}^{*})

$\theta(\mathbf{p}_n^*)$

NRH

Извините ЗНК, я неправильно понял ваш вопрос. Если вы увеличите размер выборки , смещение будет меньше, да! Оценка соответствует. Именно для равномерного распределения я бы несколько скептически фактического охвата доверительных интервалов даже при большом по причинам , которые я описал в ответе. Для всех других дистрибутивов применяется CLT, и различные методы будут производить асимптотически правильное покрытие для .

n

$n$

n

$n$

n \to \infty

$n \to \infty$

NRH

Как указывает ответ @NRH, проблема не в том, что начальная загрузка дала необъективный результат. Дело в том, что простая «подключаемая» оценка энтропии Шеннона, основанная на данных из выборки, смещена вниз от истинного значения совокупности.

Эта проблема была признана в 1950-х годах, через несколько лет после определения этого индекса. В этой статье рассматриваются основные вопросы, со ссылками на соответствующую литературу.

Проблема возникает из-за нелинейной связи индивидуальных вероятностей с этой мерой энтропии. В этом случае наблюдаемая доля генотипа для гена i в образце n , , является непредвзятой оценкой истинной вероятности, . Но когда это наблюдаемое значение применяется к формуле «подключи» для энтропии по M генам: $\hat{p}_{n,i}$ $p_{n,i}$

{\hat{θ}}_{n} = θ ({\hat{p}}_{n}) = - \sum_{i = 1}^{M} {\hat{p}}_{n, i} \log {\hat{p}}_{n, i} .

$\hat{\theta}_n = \theta(\hat{\mathbf{p}}_n) = - \sum_{i=1}^{M} \hat{p}_{n,i} \log \hat{p}_{n,i}.$

Нелинейное отношение означает, что результирующее значение представляет собой предвзятую заниженную оценку истинного генетического разнообразия.

Смещения зависит от количества генов, и число наблюдений, . Для первого порядка оценка плагина будет ниже, чем истинная энтропия на величину . Исправления более высокого порядка оцениваются в статье, приведенной выше. $M$ $N$ $(M -1)/2N$

В R есть пакеты, которые решают эту проблему. В simbootчастности, пакет имеет функцию, estShannonfкоторая вносит эти поправки смещения, и функцию sbdivдля вычисления доверительных интервалов. Лучше использовать для анализа такие устоявшиеся инструменты с открытым исходным кодом, чем пытаться начать все заново.

магистр педагогических наук
источник

Таким образом, оценка сама по себе является ошибочной из-за размера выборки? В simbootпакете выглядит многообещающе, но , кажется , не подходят для моих целей , как это нужно контрольный образец для оценки доверительных интервалов.

ZNK

«Ошибочный» не совсем верно; оценщик «смещен» в том смысле, что его ожидаемое значение не совпадает с фактическим значением населения. Это не значит, что это «ошибочно»; смещенные оценки могут быть полезны, как показано компромиссом дисперсии смещения при выборе оценок. Если simbootне отвечает вашим потребностям, Google «энтропия Шеннона смещения г» ссылки на другие пакеты R , как entropy, entropartи EntropyEstimation.

ЭдМ

Существуют дополнительные проблемы, связанные с тем фактом, что некоторые генотипы, присутствующие в популяции, вероятно, будут отсутствовать в какой-либо конкретной выборке. Некоторые из пакетов R, основанных на населении и экологии, по-видимому, имеют способы решения этой проблемы.

EdM