Подбор моделей вариограммы и ковариации

Алгоритм подбора начинается с получения предварительной оценки диапазона даннхх, и этот процесс носит название первой стадии. Используйте Zj (si )для обозначения j-го измерения переменной типа k в i-той точке пространства si.

Стадия 1

Сначала модуль Geostatistical Analyst сопоставляет каждый набор данных, Z. (si ) = Z* (si) / где s - стандартное отклонение выборки. Стадия 1 начинается с предположения, что мы имеем дело с изотропной моделью; впчисляется эмпирическая вариограмма (или ковариация) для сопоставленных даннгх Z* (si) с использованием метода сектора (как было описано ранее в разделе Бининг (объединение) оценок вариограммы и ковариации в классы лагов) для большого диапазона классов лагов, которые увеличиваются в геометрической прогрессии. Классы лагов образуются из интервалов 2, d*+1/2) , где d = 1.25, и k меняется от наименьшего аппаратного целого до наибольшего. Центр каждого класса лагов принимается равным й* cosh( logd* ) . Очевидно, что многие классы лагов останутся пустыми, и модуль Geostatistical Analyst использует только те, которые содержат данные. Назовем это эмпирической (взаимной) ковариацией Cijh*) , где i обозначает i-ый тип переменной, j обозначаетj-ый тип переменной, и k обозначает k-ый класс лагов. Первая итерация оценок параметров выполняется путем минимизации,

ZZZ-j(h*) (C~j (h; Иj)-(h ))2

i=1 j=1 к=1

для в, где в. вектор параметров для ковариационной функции i,j , и в содержит все ковариационные параметры, где

и Aij(hk) - число пар в эмпирической (взаимной) ковариационной функции для переменнгх i и j в классе k. Назовем эту оценку в(1). В следующей итерации модуль Geostatistical Analyst использует взвешенные наименьшие квадраты Кресси (1985) путем повторной минимизации (1), только в этот раз, (3)

(h к; и(;))=-л

N j (hк)

C i, (0; и(1))Cjj (0; и;) + C 2 (h *; и(;))

а затем эти веса нормируются таким образом, что каждому значению (взаимной) ковариации присваиваются равные веса,

и(1))

Назовем эту оценку в(2). Обратите внимание, что если вместо ковариации мы используем вариограмму, и(2) равно

где wi.(hk) определяется по формуле (2), а затем и2 вгчисляется по формуле (5) с весовыми значениями как в формуле (4), но теперь

j(hk;-(1))=

Оценки в(1) и в(2) - это два шага в алгоритме, который использует несколько итераций для пересчета весовых значений по методу наименьших квадратов.

Оценка в(2) используется только для оценки диапазона для размера лага, предложенного по умолчанию для метода грида, при выполнении оценки эмпирической вариограммы или ковариа-ции. Число лагов, предлагаемое по умолчанию, равно 12, следовательно, размер лага для метода грида принимается равным 2*диапазон/12 .

Шаг 2

Шаг 2 по сути повторяет шаг 1, но работает с эмпирической вариограммой или (взаимной) ковариацией сопоставленных данных Zk (si) которые используют метод грида (описанный ранее в разделе Объединение (бининг) оценок вариограммы и ковариации), где размер лага, предложенный по умолчанию, определяется по оценке диапазона из итерации в(2) первого шага. Он также допускает использование анизотропии и линейных комбинаций из не более, чем трех моделей (взаимной) ковариа-ции или вариограммы, с учетом эффекта самородка для классов лагов для каждого набора данных,

u =1

В данной формуле, Bu(i,/) - параметр частичного порога, а также i,j -тый компонент Е, а т х т положительно-определенная матрица, в которой Т - число deiia переменных, S - число различный моделей (взаимной) ковариации, используемых в линейной комбинации, а функция p (h; ф) - нормализованная модель ковариации; р (0; ф) = 1, где - параметры, которые обычно контролируют диапазон (и/или форму) модели ковари-ации. Как и ранее, в включает все параметры. Третья итерация оценок параметров, в(3) осуществляется путем минимизации формулы (1) с весами (2) для эмпирической ковариации с использованием метода грида, а затем четвертая итерация в(4) выполняется путем минимизации (1) с весами, взятыми из формул (4) и (3) для эмпирической вариограммы, использующей метод грида. Эти формулы очевидно меняются, если мы используем вариограммы, как показано для шага 1. Теперь вернемся к исходному масштабу. Окончательные модели (взаимной) ковариации выглядят следующем образом:

Cj((h) = 8. (h;и(;)),

и для вариограмм:

Если пользователь меняет какие-либо параметры (например, размер лага), оценки пересчитываются, начиная со второго шага.

Диалог двумерного распределения

Дизъюнктивный кригинг требует, чтобы все пары данных подчинялись двумерному нормальному распределению. Это условие тяжело проверить на практике. Модуль Geostatistical Analyst предоставляет визуальный инструмент, который помогает оценить выполнение условия двумерной нормальности. Теоретическая кривая, как функция лага, может быть построена на основе различных пороговых значений для функции кумулятивного распределения. (см. Deutsch и Journel, 1992, стр. 139 и Goovaerts, 1998, стр. 265). Эту теоретическую кривую можно сравнить с эмпирической вариограммой на основе индикаторов. Более широко, если графики КК показывают нормальное распределение компоненты многомерной случайной величины, и похоже, что данные подчиняются двумерному нормальному распределению, разумно предположить и наличие многомерного нормального распределения данных. Таким образом, проверка на двумерное нормальное распределение может быть использована для простого кригинга, что позволяет пользователю удостовериться в том, что карты квантилей и вероятностей основаны на разумных предположениях.

Формулы кригинга

Честь разработки методов кригинга принадлежит многим авторам. Кригинг с ковариацией эквивалентен методу лучшей линейной несмещенной интерполяции (best linear unbiased prediction - BLUP); методы простого, ординарного и универсального кригинга были разработаны такими ученхми, как Wold (1938), Холмогоров (1941), Wiener (1949), Гандин (1959), Goldberger (1962), и Henderson (1963). Методы пространственной интерполяции с использованием вариограмм были сформулированы Гандин1м (1959, 1963) и Matheron (1962, 1969). Обратитесь к работам Journel (1983) по индикаторному кригингу, Sullivan (1984) по вероятностному кригингу, и Matheron (1976) по дизъюнктивному кригингу. В работе Cressie (1990) дана подробная информация об основах кригинга.

Модуль Geostatistical Analyst включает интерполяторы, которые могут учитывать ошибку измерений. Эти модели рассматриваются в работах Гандина (1959, 1960, 1963) и Cressie (1986, 1988, 1993, стр. 127-135). К этим же моделям относятся и часто приводимые во многих учебниках по геостатистике так называемые жесткие интерполяторы. Под понятием жесткий подразумевается, что если интерполяция выполняется для точки, в которой отбирались данные, вычисленное значение будет в точности равно значению, полученному при выполнении измерений в этой точке, и для нее стандартная ошибка интерполяции равна нулю. Это может приводить к созданию странных карт, поскольку в опорных точках проинтерполированные значения будут как бы подскакивать . При наличии ошибок измерений, у пользователя может возникнуть желание отфильтровать измерения и составить более сглаженную карту проинтерполи-рованных значений. Моделирование ошибок измерений возможно только для ординарного, простого и универсального кригинга. Описание всех моделей приведено далее. Начнем с замечания. Иногда, необходимо будет иметь дело с множественными опорными точками, с множественными измерения-

ми в одной точке (ошибкой измерений) для множественных типов переменных (при кокригинге). Для обозначения j-ого измерения переменной типа k в i-ой точке пространства si мы используем обозначение (si) .

Оценка ошибки измерений

Если есть несколько наблюдений (измерений) для каждой опорной точки, модуль Geostatistical Analyst может оценить ошибку измерений. Формула для вычисления ошибки измерений следующая,

XX(Zj (si) - Z (si ))2

s, eDj=1

N - nD

где D - набор всех опорных точек, в котором выполнено более одного измерения, Zj(si.) - j-ое измерение в точке s.,

среднее значение в точке s , n

г - ME

число наблюдений в точке

D, N = Li n. для всех si в D, и nD - число опорных точек в D. Движок в диалоге по умолчанию устанавливается в положении, которое соответствует этому значению; это означает, что он установлен на значении 100 ( <С€ / самородок) процентов. Пользователи при желании могут заменить CEMe ; при выполнении оценки устанавливается, что эффект самородка не может быть меньше, чем CEME ; когда пользователь меняет значение самородка или значение, тогда сохраняется неравенство самородок > <Me