Прогнозирование экстремальных значений по выборке

[Александрович]

Имеется конечная генеральная совокупность (ГС) нормально распределенных ( с параметрами $M$ и $S$) случайных значений некоторой физической величины объемом - $N$. Из нее случайным образом извлечена выборка гораздо меньшего объема - $n$. По ней вычислены выборочные значения $m$ и $s$, которые являются эффективными, состоятельными и несмещенными оценками генеральных $M$ и $S$. Если проранжировать значения ГС в порядке возрастания, то $X_N$ - наибольшее значение в ГС. Требуется по полученным значениям $m$ и $s$ найти в аналитическом виде значение $X_{max}$ такое, чтобы $P(X_N<X_{max})=p$. Доверительная вероятность $p$ задается, пусть для определенности 0,95. Аналитический вид такой: $X_{max} = m + s*t(a,n,N)$, где $t(a,n,N)$ - некий коэффициент (не Стьюдента, но возможно связан с ним), $a = 1-p$ - уровень значимости. Трудность заключается в подборе формулы для $t(a,n,N)$, поскольку готового решения я не нашел. Прошу направить меня в правильную сторону. Спасибо.

[PAV]

PAV:

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

(формулы оформите как принято, подробности в теме)

Добавлено спустя 1 час 44 минуты 2 секунды:

Возвращено

[AndreyXYZ]

Если бы мы знали $M$ и $S$, то сразу бы могли сказать ответ.
Мы можем найди доверительные интервалы, в которых с некоторой вероятностью (которую мы зададим), будут лежать $M$ и $S$. Ход дальнейшего решения точно сказать не могу.
А откуда задачка?

[Александрович]

Интересно какой-же? Это частный случай для $n=N$. А задача практическая, из области выборочного контроля качества.

[PAV]

Если пытаться вывести разумное решение, то предлагаю начать со следующего. Предположим, что нам известны истинные значения $M$ и $S$ параметров ГС. Можете тогда выписать решение задачи через $M$ и $S$? После этого теоретически можно использовать вместо них оценки $m$ и $s$, взяв их с "достаточным запасом", определяемым подходящим доверительным интервалом.

Менее "тупо-прямолинейное" решение заключается в том, чтобы посмотреть в сторону теории порядковых статистик. Нас интересует максимум из $N$ случайных величин, что значительно проще, чем общий случай. Посмотреть, что известно про распределение. На самом деле наверняка изучались и оценки параметров, и доверительные интервалы и т.д.

Добавлено спустя 1 минуту 3 секунды:

Собственно, первая часть сообщения - это повтор того, что предложил AndreyXYZ

[Александрович]

В том-то и дело, что не могу. Через функцию распределения крайней порядковой статистистики предельная граница определяется аналитически точно только при известных параметрах математического ожидания и стандартного отклонения нормального распределения. Для этого ничего не надо измерять. А ГС это хоть и достаточно большая, но все-таки выборка с оценками параметров распределения, которому она принадлежит. Значение $t(a,n,N)$ естественно зависит от типа исходного распределения. Для определенности выбрал двухпараметрическое нормальное распределения. Если решение найдено для него, то найдено для любого вида двухпараметрических
распределений. Вопрос заключается в строгом нахождении, как Вы правильно заметили "достаточного запаса" для $m$ и $s$.

[GAA]

Мне кажется условие не очень естественное. Я предлагаю рассмотреть следующую модификацию (попутно введу обозначения с которыми, как мне кажется, будет меньше путаницы).

0. Задача. Пусть элементы «генеральной совокупности» (ГС) $X_i$, $i =1,\ldots, L+N$ — независимые, $(\mu, \sigma)$-нормально распределенные случайные величины. По выборке объемом $L$ из ГС находятся оценки $M = \frac {\sum_1^L X_i}{L}$, $S^2 = \frac {\sum_1^L (X_i - M)^2}{L-1}$. Статистика $T = \sqrt L \frac {X_{(N)} - M}{S}$, где $X_{(N)} = \max\limits_{i=L+1,\ldots, N} {X_i}$ не зависит от (неизвестных) параметров распределения, поэтому далее будет считать, что случайные величины $X_i$ имеют стандартное нормальное распределение. Для заданных $L$ и $N$ требуется найти квантили распределения $T$.
[Эта постановка задачи отличается от первоначальной тем, что «максимум берется» не по всей ГС, а по части, не попавшей в выборку. Для выборочного среднего выборки объема $L$ я использую обозначение $M$, а не $m$, а максимальную статистику обозначаю через $X_{(N)}$, а не $X_N$]

1. Решение. Обозначим функцию стандартного нормального распределения через $\Phi(z)$. Величины $X_{(N)}$, $\sqrt{L} M $, $(L-1)S^2$ независимы и, соответственно, имеют плотности $f(z) = \frac{d\Phi^N(z) }{dz}$, $f(x)=\frac{1}{\sqrt{2\pi}} \exp(-x^2/2)$ и $f(y) = \frac{1}{2^{\frac{L-1}{2}} \Gamma (\frac{L-1}{2})}y^{\frac{L-1}{2}-1}\exp(-y/2), y>0$ (см, например, [1] и [2]).
Функцию $F(t) = \mathsf P\left\{ \sqrt{L} \frac{X_{(N)} - M}{S} < t \right\}$ распределения статистики $T$ можно задать тройным интегралом

$$F(t) = \iiint\limits_{\sqrt{L}\frac{z - x}{\sqrt{y}/\sqrt{L-1}} < t}\frac{d\Phi^N(z)}{dz}\frac{1}{\sqrt{2\pi}} \exp(-x^2/2) \frac{1}{2^{\frac{L-1}{2}} \Gamma (\frac{L-1}{2})}y^{\frac{L-1}{2}-1}\exp(-y/2)dz dx dy$$

.
Сводя кратный к повторному, функцию распределения можно представить в виде

$$F(t) = \int\limits_{-\infty}^{+\infty}\frac{1}{\sqrt{2\pi}} e^{-x^2/2}dx \int\limits_0^{+\infty} \left (\int_{-\infty}^{\frac {x}{\sqrt{L}} + \frac{t\sqrt{y}}{\sqrt{L(L-1)}} }\frac{1}{\sqrt{2\pi}}e^{-u^2/2}du \right)^N \frac{1}{2^{\frac{L-1}{2}} \Gamma \left(\frac{L-1}{2}\right)}y^{\frac{L-1}{2}-1}e^{-y/2} dy$$

.
Упростить или выразить через известные функции этот интеграл я не пробовал. Квантили $t_{L, N, p}$ можно найти численно, решая уравнение $F(t_{L, N, p})= p$ для заданных $L$ и $N$, составив, тем самым, таблицу с тремя входами.
Для примера, в Maple 12 мною были найдены $t_{3, 5, 0.85} = 5.898$ и $t_{5, 10, 0.85} = 6.779$.

2. Демонстрационный пример. Пусть эксперимент включает:
0) получение выборки объемом $L+N$;
1) вычисление по подвыборке объема $L$ частных значений оценок $m$, $s$;
2) вычисление по остальной части выборки — по подвыборке объема $N$ — $x_{(N)}$ — частного значения статистики $X_{(N)}$;
3) проверку выполнения неравенства $x_{(N)} < m + \frac {t} {\sqrt L} S$ и увеличение счетчика count (числа выполнений условия), если условие выполнилось.
Тогда, при достаточно большом числе повторений (Rep) эксперимента, оценка вероятности $p^* = count/Rep$ будет приближенно равна $p$ (независимо от значений параметров распределения).
Для демонстрации этого была написан простой пример в среде TP 7.0. Для того, чтобы не загромождать текст, «генератор» нормально распределенных случайных величин вынесен в отдельный модуль NormRND. В данном примере используется две функции этого модуля: InitN1DRnd(mu, sigma) — задающая параметры распределения и инициализирующая «внутренние» переменные модуля, и N1DRnd() — возвращающая очередное значение нормально распределенной случайной величины. Модуль реализует, изложенный в [3], метод «полярных координат». В качестве базового генератора равномерного распределения используется генератор модуля system.

uses NormRND;
const
 p = 0.85;
 t : Double = 5.8976;
 L = 3;
 N = 5;
 mu: Double = 1;
 sigma: Double  = 1;
 Rep = 1000000;
var
 X: array[1..L+N] of Double;
 m, s, xN: Double;
 Count: LongInt;
 i, j: LongInt;
 C:    Double;
begin
 C:= t/sqrt(L);
 Count:= 0;
 InitN1DRnd(mu, sigma);
 for j:= 1 to Rep
  do begin
      for i:= 1 to L+N do X[i]:= N1DRnd;
      m:= 0; for i:= 1 to L do m:= m + X[i]; m:= m/L;
      s:= 0; for i:= 1 to L do s:= s + sqr(X[i]-m); s:= sqrt(s/(L-1));
      xN:= X[L+1]; for i:= L+2 to L+N do if X[i] > xN then xN:= X[i];
      if xN < m + C*s then inc(Count);
     end;
 writeln('p* =',  Count / Rep);
end.

Для $L = 3$, $N = 5$, $Rep = 1000000$ и $p = 0.85$ получено $p^* \approx 0.85$ при различных значениях $\mu$ и $\sigma$. Аналогичный результат был получен для $L = 5$, $N = 10$.

ref
[1] Приложение «Минимум и максимум» к лекциям Н.И. Черновой по ТВ за 2005 г. Главу 6 «Распределения связанные с нормальным» лекций Н.И. Черновой по МС за 2006 г.
[2] Крамер Г. Математические методы статистики. — М.: Мир. 1975. [По поводу независимости M и S^2 можно посмотреть]
[3] Кнут Д. Искусство программирования. т.2 (Получисленные алгоритмы).

25.01.09 исправлена опечатка в тройном интеграле: Ф(x) заменено на Ф(z)
26.01.09 исправлены еще четыре очевидные опечатки и ошибки форматирования

[Александрович]

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

[Александрович]

Квантили $t_{L,N,p}$ получились огромными, возможно из-за малости $L$. Вообще-то в практической задаче $L=30$ а $N=270$. Если Вам не трудно, уважаемый $GAA$, решите эту задачу для указанных значений $L,N$ при $p=0,75$. Вы взяли произведение трех плотностей, это - плотность распределения вероятности какого события? Вы где-нибудь встречали решение этой задачи? Если нет, то опубликуйте результаты, что бы можно было на Вас сослаться. Это актуально в выборочном разрушающем контроле качества. Хорошо бы в ГОСТ протащить с готовой простой формулой для $t_{L,N,p}$, либо в виде таблицы. И наконец, если Вы примените такой-же подход к выборочным среднему и дисперсии для нахождение доверительного интервала математического ожидания, у Вас должно получиться распределение Стьюдента. Попробуйте, это гарантия правильной идеологии при решении задачи о максимальном значении.

[GAA]

С увеличением $L$ и $N$ объем вычислений значительно увеличивается; увы, времени на выполнение Вашего пожелания у меня нет. (Я специально подробно выписал интеграл, чтобы желающие получить значения квантилей смогли это проделать самостоятельно.)
Записанное произведение трех плотностей — это плотность совместного распределения величин $X_{(N)}$, $\sqrt{L}M$, $(L-1)S^2$ (с учетом того, что $\mu=0$, $\sigma=1$).
Точно то, что я изложил в своём предыдущем сообщении, я не встречал, и публиковать пока тут нечего: написанное мгновенно следует из стандартного вузовского курса ТВиМС. Создание таблицы и/или подбор приближенного выражения для больших $N$ — дело довольно хлопотное, браться за это у меня сейчас времени нет.
«Применить такой же подход к выборочному среднему и дисперсии для нахождения доверительного интервала для математического ожидания» не смогу — его уже применили задолго до меня, и я, будучи студентом, а затем преподавателем, неоднократно воспроизводил «применение этого подхода».

[Александрович]

Спасибо, а почему

GAA писал(а):

$\sqrt{L}M$

, а не просто $M$. И все-таки куда можно сослаться, для указания наиболее близкого или подобного решения. А публиковать есть что, во-первых это новый метод обнаружения выбросов, во-вторых построение доверительного интервала для случайной величины (границы прогнозной модели или карты качества).

[Александрович]

Я понял почему получились большими $t(L,N,p)$. Вы считаете

GAA писал(а):

Статистика $T = \sqrt L \frac {X_{(N)} - M}{S}$

, то есть используете среднеквадратичное отклонение для среднего, а естественней использовать среднеквадратичное для случайной величины, то есть без умножения на \sqrt L \.

[GAA]

1. Я использовал статистику $\sqrt{L}M$, поскольку она имеет стандартное нормальное распределение (при сделанных предположениях о $\mu$ и $\sigma$), а также потому, что «Стьюдентовское отношение» имеет вид $t= \sqrt{L} \frac{M-\mu}{S}$. Думал, все это облегчит понимание того, как получено выражение для функции распределения $T$.

2. Повторюсь. Не надо ни на что ссылаться «для указания наиболее близкого или подобного решения». Выписанное мною выражение для функции распределения очевидено. Просто откройте учебники по МС — два таких учебника ([1], [2]) я привел выше, — и разберитесь с обязательным материалом. В научных публикациях не принято давать ссылки на общеизвестные факты.

3. Замечу, что доверительные области строятся для детерминированной величины. Обычно рассматривают доверительные области для параметрической (возможно векторной) функции, в частности, рассматривают доверительные интервалы для скалярного параметра. Для случайной величины доверительный интервал построить нельзя — нельзя просто по определению доверительного интервала.

4. Выбросами никогда не занимался. Но «стьюдентизированное» экстремальное отклонение устроено совсем по другому:
$T_n = \frac {X_{(n)} - M_n}{S_n}$,
где $X_{(n)}$, $M_n$ и $S_n$ вычисляются по одной и той же выборке объема $n$. (См. [4], где приведены ссылки на старинные работы, посвященные критериям исключения грубых наблюдений. См. также [5, n. 4.2.4. Нормальные экстремальные отклонения (с. 173)], [6, Выпадающие наблюдения, с. 707].)

Ref. (ссылки [1]—[3] даны ранее в этой темы)
[4] Большев Л.Н., Смирнов Н.В. Таблицы математической статистики. — М.: Наука, 1983; Критерии исключения резко выделяющихся наблюдений (c.62).
[5] Гумбель Э. Статистика экстремальных значений. — М.: Мир, 1965.
[6] Кендалл М., Стьюарт А. Статистические выводы и связи. — М.: Наука, 1973.

[Александрович]

GAA в сообщении #182957 писал(а):

3. Замечу, что доверительные области строятся для детерминированной величины. Обычно рассматривают доверительные области для параметрической (возможно векторной) функции, в частности, рассматривают доверительные интервалы для скалярного параметра. Для случайной величины доверительный интервал построить нельзя — нельзя просто по определению доверительного интервала.

Согласен, но ведь Вы показали как рассчитать $ 0,85 $ доверительную границу интервала случайной величины - максимального значения выборки. Квантиль, величина уже не случайная.

[GAA]

Ни о какой доверительной границе я не писал. Еще раз: нельзя построить доверительный интервал для случайной величины. Разберитесь, пожалуйста, с определением доверительного интервала.