Опыт применения геостатистического анализа при исследовании ресурсных зон поселений эпохи раннего средневековья в Кисловодской котловине

Д.С. Коробов

Институт археологии РАН, Москва

Работа выполнена в рамках проекта РФФИ № 18-09-00615а.

Исследования ресурсных зон поселений разных эпох и культур (off-site archaeology) становятся достаточно широко применяемой практикой в зарубежной (преимущественно британской и нидерландской) археологии, начиная с 1980-х гг. Работами Т. Уилкинсона, Т. Уильямсона, Дж. Бинтлиффа [Williamson, 1984; Wilkinson, 1989; Bintliff, 2000] и др. выработаны методы оценки использования ресурсных зон вокруг поселений, которые базируются на разработках палеоэкологической школы Э. Хиггза [Jarman, 1972; Jarman et al., 1972; Higgs, Jarman, 1975; Barker, 1975; Hodder, Orton, 1976; Foley, 1977], восходящей к концу 1960-х гг. (Site Catchment Analysis). Согласно данной теории, земледельческое население использовало окружающие ресурсы вокруг рядовых сельских поселений на удалении не более 5 км, что соответствует расстоянию в 1 час пешей ходьбы по непересеченной местности. При этом зона земледельческого использования, как правило, лежит на удалении не более 1 км от поселения.

Проверка данных теоретических положений осуществлялась зарубежными коллегами в основном методом систематического сбора подъемного материала и изучением содержания фосфатов в почвах вокруг поселений [Sjöberg, 1976; Widgren, 1983; Thurston, 2001]. Широко используемые в зарубежной археологии подходы и методы изучения ресурсных зон были адаптированы и развиты в ходе многолетнего изучения следов земледелия в Кисловодской котловине методами археологического почвоведения, итогом чему стал ряд статей и монографий [Борисов, Коробов, 2013; Korobov, Borisov, 2013]. Коллективом почвоведов из Института физико-химических и биологических проблем почвоведения РАН (г. Пущино) под руководством к.б.н. А.В. Борисова был разработан ряд уникальных индикаторов (уреазная активность почв, количество термофильных микроорганизмов), позволяющих выявить присутствие органических удобрений в почве пахотных горизонтов, описание которых также было опубликовано [Чернышева и др., 2016; Chernysheva et al., 2015; 2017]. Полученные результаты легли в основу авторской модели ресурсных зон вокруг раннесредневековых поселений Кисловодской котловины, моделирование которых осуществлялось с помощью пространственного ГИС-анализа [Коробов, 2017а; 2017б].

Настоящая работа посвящена продолжению этих исследований, выполненному на более детальном уровне. Новые результаты по моделированию активных зон древнего и средневекового земледелия были получены в ходе систематического обследования окрестностей укрепления Подкумское 3, относящегося к эпохе раннего средневековья.

1. Систематическое обследование окрестностей укрепления Подкумское 3. Памятник располагается на левом берегу р. Подкумок, в 2 км к юго-западу от окраины сел. Терезе. Он был обследован нами дважды, в 2009 и 2012 гг., тогда же проводились его шурфовка и топографическая съемка [Коробов, 2017в. С. 60–66]. В окрестностях крепости закладывалась небольшая серия почвенных разрезов на разном удалении от памятника, информация по которым послужила основанием реконструкции границ земледельческой зоны аланского поселения V–VIII вв. [Борисов, Коробов, 2013. С. 156–159; Чернышева и др., 2016. С. 102–155; Chernysheva et al., 2017].

Почвенно-археологические исследования сезона 2018 г. заключались в систематическом обследовании всей прилегающей к поселению территории путем заложения серии почвенных разрезов по сетке с шагом в 50 м. Поскольку территория, охваченная почвенно-археологическими исследованиями, по площади превышала 1 кв. км, были использованы современные возможности низковысотной аэрофотосъемки для получения трехмерного изображения поверхности с помощью фотограмметрии. Работы по созданию микротопографической модели поверхности проводились Ю.М. Свойским и Е.В. Романенко (Институт классического Востока и античности ВШЭ, Лаборатория дистанционного зондирования и анализа пространственных данных). Была получена цифровая модель поверхности с пространственным разрешением 5 см, которая была соотнесена с мировыми географическими координатами с помощью GNSS-приемника (рис. 1).

Всего в окрестностях укрепления Подкумское 3 было заложено 99 почвенных разрезов по координатной сети, которая получила буквенно-цифровое обозначение (линии разрезов в направлении юго-восток – северо-запад получили наименования латинскими буквами от A до P, по линии северо-восток – юго-запад – цифрами от 1 до 28). Начало сети (разрез А-1) устроено в 50 м к юго-западу от крайних сооружений укрепления, видимых на поверхности. При этом почвенные разрезы закладывались с максимальной плотностью (через 50 м) на прилегающей к поселению территории, где не было распашки в советское время. Далее, на удалении 400–600 м от укрепления на залежном участке разрезы закладывались с шагом в 100 м, а на удалении свыше 750–1000 м – с шагом в 200 и 400 м. Таким образом, крайняя линия разрезов B-28 – L-28 находилась на расстоянии 1400–1500 м от укрепления Подкумское 3 (рис. 1).

Почвенные разрезы размерами 1 × 1 м копались пластами по 10 см с отбором археологического материала из каждого пласта. Кроме фиксации археологических находок, из каждого пласта отбиралась серия почвенных проб для определения содержания фосфатов, магнитной восприимчивости, активности фермента уреазы, обилия кератинолитических грибов и численности термофильных бактерий. В настоящий момент перед коллективом стоит задача картографирования этой обширной информации, для чего предполагается использовать возможности геостатистического анализа, которые предоставляются современными ГИС. Целью настоящей публикации является рассмотрение основных методов подобного анализа и их адаптации для использования при изучении ресурсных зон поселений в окрестностях Кисловодска.

2. Основы геостатистического анализа и алгоритмы пространственной интерполяции. Одной из наиболее широко распространенных задач, встающих перед исследователем при изучении пространственных особенностей того или иного явления, является построение непрерывной поверхности на основе точечной информации. Простейшим примером подобного моделирования может быть создание цифровой модели рельефа местности по измерениям высоты, которые производятся с помощью геодезических приборов. Однако пространственное моделирование не ограничивается созданием математического представления о рельефе – подобное выражение может получить любой количественный признак, измеренный в ходе наблюдений. Создание непрерывной поверхности и моделирование пространственного распределения количественной информации могут осуществляться с помощью ряда методов геостатистического анализа. В настоящее время геостатистика успешно используется в географии, геофизике, изучении окружающей среды и экологии, геологии и геологоразведке, эпидемиологии и социологии, рыбном и лесном хозяйстве и т.п. [Демьянов, Савельева, 2010. С. 8]. Есть примеры использования геостатистических процедур и в археологии [Lloyd, Atkinson, 2004]. В отечественной науке известны примеры адаптации методов геногеографии для картографирования частотных характеристик археологических объектов и артефактов [Балановская, Балановский, 2007. С. 247–252; Ковалевская, 2016].

Широкий инструментарий по использованию геостатистического анализа предлагается пользователям пакета программ ArcGIS, в котором подобные процедуры пространственного моделирования осуществляются с помощью специального модуля геостатистического анализа (Geostatistical Analyst). Ниже рассматривается серия карт, построенных с помощью данного модуля при использовании лицензионной версии программы ArcGIS 10.5.1.

Адекватность полученного результата геостатистического моделирования исходным данным зависит от правильности выбора метода интерполяции имеющейся информации на все анализируемое пространство. Согласно общепринятому определению, приводимому Википедией, интерполяцией называют способ нахождения промежуточных значений величины по имеющемуся дискретному набору известных значений. Модуль Geostatistical Analyst программы ArcGIS предлагает широкий спектр подобных методов, описание которых можно найти в руководстве пользователя данным дополнительным модулем на сайте ESRI. Более подробную информацию о некоторых методах интерполяции можно получить в специальной литературе, например в монографии [Демьянов, Савельева, 2010] или широко известных учебниках по применению ГИС в археологии [Wheatley, Gillings, 2002. P. 174–177; Conolly, Lake, 2006. P. 90–102]. Близкие приводимым ниже примерам сравнения методов интерполяции при картировании мощности гумуса и количества фосфатов в почве приводятся в статьях [Мыслыва и др., 2017а; 2017б]. Следует подчеркнуть, что мною не ставилась цель математического обоснования выбора того или иного алгоритма интерполяции. Сравнение интерполяционных методов ведется на самом общем уровне, без углубления в достаточно сложные для понимания основы геостатистического моделирования, познакомиться с которыми читатель может при обращении к цитируемой выше литературе. При описании методов интерполяции использовалось вышеупомянутое руководство пользователя модулем Geostatistical Analyst.

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

Детерминированные методы интерполяции можно разделить на две группы: глобальные и локальные. Глобальные методы вычисляют проинтерполированные значения на основании всего набора данных. Локальные методы вычисляют проинтерполированные значения на основании измеренных точек в пределах окрестностей, которые являются меньшими пространственными областями внутри большей изучаемой территории. Geostatistical Analyst предоставляет метод глобального полинома (англ. Global Polynomial) в качестве глобального интерполятора, а также методы обратного взвешенного расстояния (ОВР, англ. Inverse Distance Weighting, IDW), локального полинома (англ. Local Polynomial), радиальной базисной функции (РБФ, англ. Radial Basis Functions, RBF), сглаживания ядра (англ. Kernel Smoothing) и ядра диффузии (англ. Diffusion Kernel) в качестве локальных интерполяторов.

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

Приведем примеры использования методов детерминированной интерполяции. В качестве исходных данных для моделирования используется информация о глубине почвенных разрезов, устроенных в окрестностях укрепления Подкумское 3. Данный показатель варьирует от 10 до 140 см (рис. 2). Таким образом, ниже приводится серия карт моделирования мощности почвенного покрова, построенная с помощью разных методов интерполяции.

Интерполяция по методу глобального полинома подбирает к входным опорным точкам сглаженную поверхность, определяемую математической функцией (полиномом). Поверхность глобального полинома постепенно изменяется и характеризует грубую структуру в данных. Можно задавать количество изменений поверхности вручную от 1 до 10. В результате интерполяции по методу глобального полинома формируется сглаженная поверхность, которая отражает постепенные тренды поверхности в исследуемой области (рис. 3).

Интерполяция по методу глобального полинома используется в следующих случаях:

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

2) Исследование и/или удаление эффектов трендов с большим диапазоном или глобальных трендов. При таких обстоятельствах метод часто называют анализом поверхности тренда.

Интерполяция по методу глобального полинома создает медленно изменяющуюся поверхность с помощью полиномов низкого порядка, которые, возможно, описывают некий физический процесс, например, загрязнение и направление ветра. Однако следует отметить, что чем более сложный используется полином, тем труднее приписать ему физический смысл (рис. 4). Кроме того, вычисленные поверхности весьма чувствительны к выпадающим значениям (чрезвычайно высоким и/или низким) (рис. 5).

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

Интерполяция по методу локальных полиномов будет более точной, если данные обладают следующими свойствами:

1) Образцы были отобраны по регулярной сетке, то есть через равные расстояния.

2) Значения данных в окрестности поиска нормально распределены.

На рис. 6 приведен пример моделирования почв по глубине с использованием метода локального полинома со значением исследовательского анализа поверхности тренда 44, который задается в окне использования программы и может варьировать от 1 до 100. Данный инструмент позволяет управлять одновременно всей совокупностью параметров, необходимых для анализа.

Метод обратных взвешенных расстояний (ОВР) однозначно предполагает, что объекты, которые находятся поблизости, более подобны друг другу, чем объекты, удаленные друг от друга. Чтобы проинтерполировать значение для неизмеренного положения, ОВР использует измеренные значения вокруг интерполируемого местоположения. Наиболее близкие к проинтерполированному местоположению измеренные значения оказывают большее влияние на прогнозируемое значение, чем удаленные от него на значительное расстояние. Метод ОВР предполагает, что каждая измеренная точка оказывает локальное влияние, которое уменьшается с увеличением расстояния. Это придает больший вес точкам, расположенным ближе всего к интерполируемому местоположению. Вес точки уменьшается как функция от расстояния. Поэтому метод носит название обратных взвешенных расстояний.

Наилучшие результаты при использовании метода ОВР могут быть получены в том случае, если сеть опорных точек достаточно плотная для того, чтобы отразить локальную вариацию, которую требуется смоделировать. Если сеть входных опорных точек редкая или нерегулярная, результирующая поверхность будет недостаточно адекватно представлять требуемую поверхность. Это очевидно из рис. 7, на котором представлены результаты моделирования мощности почвенного слоя возле укрепления Подкумское 3. Следует отметить, что результат построения поверхности может зависеть от величины радиуса окрестностей поиска интерполируемых значений. В нашем случае здесь и ниже данный радиус был равен 50 м – шагу систематической сетки почвенных разрезов, устроенных в ходе полевых работ.

Радиальные базисные функции (РБФ) представляют собой еще один набор методов жесткой интерполяции; это означает, что поверхность должна проходить через каждое измеренное опорное значение. Существует пять различных базисных функций: плоский сплайн, сплайн с натяжением, полностью регуляризованный сплайн, функция мультиквадриков и функция обратных мультиквадриков. Каждая базисная функция имеет различную форму и позволяет получать разные интерполированные поверхности. В нашем случае использовалась функция сглаживания поверхности в виде полностью регуляризованного сплайна (рис. 8).

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

Радиальные базисные функции используют для создания сглаженных поверхностей из большого количества расчетных данных. С помощью функции можно успешно создавать слабо изменяющиеся поверхности, например поверхности высот. Однако эти методы не подходят при значительных изменениях значений поверхности в пределах коротких расстояний или в том случае, если опорные данные могут быть ошибочными или неопределенными. В нашем случае модели мощности почвенного слоя, полученные с помощью методов ОВР и РБФ, оказались весьма близкими (рис. 7, 8).

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

Интерполяция диффузии основывается на фундаментальном решении уравнения переноса тепла, которое описывает, как тепло или частицы со временем распределяются в однородной среде. Проинтерполированные значения, полученные с помощью этого метода, плавно «обтекают» барьеры. В отсутствие барьеров интерполированные значения, полученные по методу интерполяции диффузии, приблизительно соответствуют интерполированным значениям, полученным методом интерполяции гауссова ядра (рис. 9). В программе имеется возможность задать ширину анализируемой полосы, в зависимости от которой меняются полученные результаты пространственного моделирования (рис. 10).

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

Интерполяция ядра использует следующие радиально-симметричные ядра: экспоненциальное, гауссово, четвертого порядка, Епанечникова, полиномиальное пятого порядка и константу. Ширина ядра определяется прямоугольником вокруг наблюдений. В качестве примеров приведены результаты моделирования мощности почвенного слоя, полученные с использованием интерполяции экспоненциального (рис. 11) и гауссова ядра (рис. 12), результаты которых выглядят практически идентичными.

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

В геостатистике предполагается, что изучаемая территория является результатом случайного процесса. Геостатистика выполняет две ключевые задачи: раскрытие правил зависимостей и выполнение интерполяции. Прежде чем приступать к интерполяции, необходимо изучить правила зависимостей.

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

Многие методы связаны с геостатистикой, но все они принадлежат семейству методов кригинга. Данный интерполяционный метод геостатистики назван в честь южноафриканского горного инженера Дэниела Крига, занимавшегося ручным созданием геологических карт по ограниченному набору данных в изучаемой области.

При использовании модуля ArcGIS Geostatistical Analyst доступны следующие методы кригинга: ординарный, простой, универсальный, вероятностный, индикаторный и дизъюнктивный (англ. Ordinary, Simple, Universal, Probability, Indicator, Disjunctive), а также их аналоги в кокригинге. Эти методы кригинга не только создают поверхности интерполяции и ошибок, но также строят выходные карты вероятности и квантиля, в зависимости от пользовательских потребностей.
Кригинг имеет две основные цели: количественной оценки пространственной структуры данных и создания интерполяций. При количественной оценке структуры, которая еще называется вариографией, подбирается модель пространственной зависимости к своим данным. Чтобы выполнить интерполяцию неизвестного значения для определенного местоположения, метод кригинга использует подобранную в вариографии модель, конфигурацию пространственных данных и значения измеренных опорных точек вокруг прогнозируемого местоположения. Модуль ArcGIS Geostatistical Analyst оснащен инструментами, которые помогают определять необходимые параметры, а также предоставляет надежные значения по умолчанию, с помощью которых можно быстро создать поверхность.
Кригинг основан на двух задачах: во-первых, оценка функций вариограммы и ковариации для значений статистической зависимости (называемой пространственной автокорреляцией) и, во-вторых, интерполяция неизвестных значений с использованием обобщенных линейных методов регрессии (кригинга). Это две разные задачи, поэтому принято считать, что в геостатистике данные используются дважды: вначале для оценки пространственной автокорреляции, а затем для выполнения интерполяции.

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

Ординарный (обычный) кригинг – более широко используемый из методов кригинга. Основой данного метода является предположение, что постоянное среднее значение неизвестно. В ординарном кригинге учитываются не только расстояния от интерполируемой точки, но и расстояние между самими точками так, что вес более близких друг к другу точек уменьшается (рис. 13).

В отличие от ординарного, простой кригинг предполагает, что постоянное среднее значение известно. Использование в качестве среднего его статистической оценки (математического ожидания) делает веса зависимыми от значений исходного набора данных (рис. 14). Кроме того, оценка математического ожидания может оказаться искаженной, смещенной и т.п. Поэтому простой кригинг редко применяется как самостоятельный метод оценивания, обычно его использование связано с искусственными комбинациями, где среднее известно вследствие предварительных манипуляций с исходными данными [Демьянов, Савельева, 2010. С. 115–116].

При использовании метода универсального кригинга предполагается, что в данных имеется какая-либо доминирующая тенденция, которую можно смоделировать с помощью детерминистической полиномиальной функции. Данный полином вычисляется из исходных значений измерений, и автокорреляция моделируется по случайным ошибкам. Когда к случайным ошибкам подобрана модель, перед расчетом полином суммируется с полученными данными, чтобы получился осмысленный результат. Универсальный кригинг применяется тогда, когда известно, что в данных существуют определенные тенденции (рис. 15) и имеется возможность провести научное описание для их подтверждения.

Метод индикаторного кригинга использует модель, по которой непрерывные данные преобразуются в двоичные переменные с использованием вычисленного порогового значения. Построенные модели будут выделять области превышения данного порогового значения (рис. 16) и области, где пороговое значение не превышено (рис. 17). В нашем случае пороговым значением стала мощность почвенного покрова около 1 м.

Существует ряд других методов кригинга, основанных на применении вероятностных процедур (вероятностный кригинг, дизъюнктивный кригинг – рис. 18), использование которых требует специальной подготовки. Для облегчения использования метода кригинга и автоматизации сложного процесса геостатистического моделирования в модуле Geostatistical Analyst разработана процедура эмпирического байесовского кригинга.

Эмпирический байесовский кригинг (ЭБК, англ. Empirical Bayesian kriging, EBK) – это метод геостатистической интерполяции, автоматизирующий наиболее трудоемкие аспекты построения корректной модели кригинга. Другие методы кригинга в модуле Geostatistical Analyst требуют ручного изменения параметров для достижения точных результатов, в то время как метод ЭБК автоматически вычисляет эти параметры путем разбиения данных на поднаборы и их моделирования (рис. 19).

Кроме того, эмпирический байесовский кригинг отличается от других методов кригинга тем, что учитывает ошибку, связанную с оценкой основной вариограммы. Другие методы кригинга рассчитывают вариограмму на основе известных местоположений данных и используют эту единственную вариограмму для прогнозирования в неизвестных местоположениях; данный процесс неявно предполагает, что расчетная вариограмма является истинной для региона интерполяции. Не учитывая неопределенность расчета вариограммы, другие методы кригинга недооценивают стандартные ошибки интерполяции.

Полученные в ходе сравнительного анализа модели мощности почвенного покрова в окрестностях укрепления Подкумское 3 показывают близкие результаты при использовании алгоритмов ОВР и РБФ, а также простого, дизъюнктивного и эмпирического байесовского кригинга (рис. 7, 8, 14, 18, 19). Учитывая сложность математической процедуры построения геостатистических моделей интерполяции, было решено остановиться на методе обратных взвешенных расстояний как наиболее пригодном для работы с регулярной сеткой измеренных данных. Полученные ниже карты некоторых вычисленных характеристик ресурсной зоны окрестностей укрепления Подкумское 3 были построены с использованием данного алгоритма.

3. Анализ распределения керамики в окрестностях укрепления Подкумское 3. Керамические фрагменты, обнаруживаемые вне поселений, прекрасно маркируют зону сельскохозяйственного освоения (off-site), поскольку являются хорошим индикатором внесения органических удобрений на поля, в процессе которого вместе с навозом в пахотные слои поступал бытовой мусор [Williamson, 1984; Wilkinson, 1989; Miller, Gleason, 1994. P. 37–38; Ford et al., 1994; O’Connor, Evans, 2005. P. 245]. Помимо керамики, в качестве индикаторов внесения удобрений могут выступать повышенное содержание фосфатов, высокие показатели уреазной активности и значительное количество термофильных бактерий [Чернышева и др., 2016. С. 110–155; Chernysheva et al., 2015; 2017]. Эти показатели будут непременно учтены в ходе дальнейшей работы над моделированием ресурсной зоны укрепления Подкумское 3. В настоящей статье речь пойдет только о находках керамики как маркере сельскохозяйственной активности древнего и средневекового населения данного памятника.

Практически во всех разрезах присутствовала керамика, в некоторых из них – кости животных. Всего было обнаружено около 2380 фрагментов сосудов, из которых около 800 фрагментов были атрибутированы как относящиеся к аланской культуре эпохи раннего средневековья, более 1300 – к кобанской культуре позднего бронзового – раннего железного века; около 270 фрагментов остались неопределенными, но с высокой долей вероятности датирующимися первыми веками н.э.

При этом степень фрагментированности керамики разных эпох и культур может весьма сильно различаться. Так, в среднем черепки раннесредневековой керамики приблизительно в два раза крупнее, чем фрагменты сосудов кобанской культуры. Для достижения большей объективности при геостатистическом сопоставлении пространственного распределения керамики было решено использовать не количество обнаруженных фрагментов, а их вес. Общий вес керамических фрагментов, найденных в разрезах 2018 г. в окрестностях укрепления Подкумское 3, составил около 9 кг. Керамика кобанской культуры, составлявшая более половины от обнаруженной по количеству, весит чуть более 4 кг. Аланская керамика (примерно треть от найденной в 2018 г.) имеет общий вес около 3,6 кг, а неопределенная по культурной принадлежности керамика – чуть более 1,1 кг.

Если обратиться к общему распределению керамики по разрезам, то становится очевидным ее неравномерное залегание на обследованной территории. Так, керамика кобанской культуры имеет максимальные значения веса в двух разрезах – А-8 и D-28, где было найдено более 500–600 гр. посуды этого периода (рис. 20). Высокие значения веса кобанской керамики в последнем разрезе в сочетании со значительным количеством ее (около 100 гр. и выше) в разрезах B-20, B-28 и D-20 создают обширный ареал следов кобанского земледелия в южной части обследуемой территории. Следует отметить, что расстояние между рядами разрезов здесь составляет 400 м, и смоделированная обширная территория высоких значений веса керамики кобанской культуры является результатом интерполяции.

Гораздо более обоснованным выглядит ареал сельскохозяйственной активности населения кобанской культуры на территории, примыкающей к укреплению Подкумское 3 (линии разрезов 4–12) (рис. 20). Обращает на себя внимание его ограниченность с востока, севера и запада. Максимальные значения веса кобанской керамики приходятся на разрез А-8, где в нижних пластах был обнаружен развал крупного сосуда. Очевидно, в данном случае мы имеем дело со следами небольшого поселения кобанской эпохи, которое прекрасно оконтуривается по результатам геостатистического анализа (рис. 20). К поселению примыкает зона земледельческой активности, следы которой прослеживаются в разрезах в виде находок кобанской керамики на расстоянии 300–350 м. Отдельный «выплеск» с относительно высокими значениями веса кобанской керамики (74 гр.) наблюдается в разрезе J-8, находящемся на расстоянии около 450 м от разреза А-8 с предполагаемыми слоями поселения. Далее значения веса керамики снижаются до 10–30 гр. на разрез, что очевидно маркирует границу ресурсной зоны означенного кобанского поселения. Большое количество находок кобанской керамики в крайних разрезах линии 28, вероятно, относится уже к ресурсной зоне другого поселения (рис. 20), земледельческие угодья которого сопровождались участками террасного земледелия (на них приходятся разрезы B-28 – F-28) классического облика для памятников кобанской культуры Кисловодской котловины IX–VI вв. до н.э. [Борисов, Коробов, 2013. С. 68–98].

Иной ареал сельскохозяйственной активности прослеживается при картографировании фрагментов посуды эпохи раннего средневековья (рис. 21). Здесь максимальные значения веса керамики наблюдаются в разрезах, непосредственно примыкающих к укреплению Подкумское 3 (А-1 – А-4) с максимальным количеством ее в разрезе А-3 (655 гр.). Этот разрез явно маркирует зону аланского неукрепленного поселения, расположенного за пределами видимой на поверхности укрепленной части. Далее следует ресурсная зона интенсивной сельскохозяйственной деятельности, которая выявляется за счет находок аланской керамики в разрезах, где ее вес колеблется в пределах 50–100 гр. на разрез (рис. 21). Границы этого ареала проходят на расстоянии 500–600 м от поселения, что прекрасно соотносится с теорией ресурсных зон Э. Хиггза и его школы [Jarman, 1972; Jarman et al., 1972; Higgs, Jarman, 1975; Barker, 1975; Hodder, Orton, 1976. P. 229–236; Foley, 1977]. Интересно отметить, что за понижением количества найденной аланской керамики до минимальных значений около 10 гр. следует участок с высокими показателями веса керамических фрагментов (60–80 гр.), маркируемый разрезами J-16 – N-16, который отстоит от поселения на расстояние около 1000 м. Вероятно, здесь мы имеем дело с регулярно обрабатывавшейся заимкой, также относящейся к укрепленному поселению Подкумское 3, пространство между которой и основным ареалом сельхозугодий находилось в менее интенсивной обработке (рис. 21).

Наблюдающееся в юго-западной части исследованной разрезами территории повышенное содержание аланской керамики (разрезы D-28 – H-28, 50–70 гр. на разрез), по аналогии с описанной выше кобанской, очевидно относится уже к ресурсной зоне другого поселения, расположенного далее к юго-западу.

Сложнее для интерпретации распределение керамики, определение которой затруднительно. Очевидно, сама процедура отнесения того или иного черепка к данной категории находок, описанная ранее [Коробов, 2017б. Т. 1. С. 167], приводит к тому, что эта посуда занимает промежуточное положение между керамикой кобанской и аланской культур, имея схожие черты и с той, и с другой. Косвенно это подтверждается результатами геостатистического анализа неопределенной керамики, происходящей из 59 почвенных разрезов (рис. 22). Наибольшее ее количество приходится на разрез B-10, где обнаружено 166 гр. подобной керамики. Более 50 гр. найдено в разрезе D-10, находящемся неподалеку. В остальных случаях в разрезах встречается менее 50 гр. такой керамики, причем ареал ее совпадает с одной стороны с ареалом находок кобанской посуды – он отстоит от укрепленного поселения и занимает центральную часть мыса между линиями разрезов 5 и 12. Некоторое количество неопределенной керамики встречается в разрезе B-28 (33 гр.), где кобанская посуда встречается в изобилии. С другой стороны, в разрезах J-16 и L-16, где обнаружено значительное количество аланской керамики при практически полном отсутствии кобанской, имеется и 30–50 гр. керамики неопределенной. Таким образом, не исключено, что схожесть неопределенной керамики и с кобанской, и с аланской посудой влияет на полученный результат геостатистического анализа ее пространственного распределения.

Интересно проследить залегание посуды разных культурно-хронологических этапов по глубине. Поскольку керамика в разрезах отбиралась с каждых 10 см, имеется возможность отобразить эту информацию картографически. Попластовое распределение кобанской керамики (анимационный ролик 1) демонстрирует ее незначительное количество в верхних двух пластах; на третьем и четвертом пластах максимальное количество керамики зафиксировано в разрезах B-20 и D-20. Далее на большей глубине (пласты 5–8) самое большое количество керамики – от 55 до 185 гр. – найдено в разрезе D-28, находящемся на древней террасе. Здесь кобанская керамика залегает в слое погребенной почвы. На восьмом пласте в разрезе А-8 обнаружен развал крупного сосуда весом в 385 гр., что также отразилось на результатах геостатистического анализа. Наконец, девятый пласт дает значительное количество керамики в разрезе B-6 (75 гр.) и небольшое количество посуды в нижних слоях погребенной почвы на террасах в разрезах B-28 и D-28 (23–26 гр.).

Иные закономерности прослеживаются при попластовом анализе распределения аланской керамики (анимационный ролик 2). Здесь максимум керамики встречается в верхних трех пластах, в особенности в примыкающих к укреплению разрезах A-1 – A-3. Второй и третий пласты, очевидно, маркируют залегание аланского пахотного слоя, поскольку именно в них повсеместно встречается керамика раннего средневековья. В некоторых случаях этот слой, скорее всего, продолжается и на четвертом пласте. По крайней мере, в разрезах J-8, N-16 и J-16 именно на этой глубине встречено от 32 до 64 гр. аланской керамики. Как уже упоминалось выше, значительное количество керамики из разреза А-3 (более 100 гр.) относится, по-видимому, к культурному слою существовавшего здесь поселения. Далее на пятом пласте количество керамики резко снижается – она встречена в значительном количестве лишь в разрезе А-4 (более 100 гр.), который также, скорее всего, относится к поселению. Ниже количество аланской керамики минимально (немного ее найдено лишь в разрезе В-28 – 22 гр.), а начиная с седьмого пласта, она исчезает вовсе.

 Как и в описанном выше случае с неопределенной керамикой, ее распределение по пластам напоминает ситуацию и с кобанской, и с аланской посудой. Максимальное количество неопределенной керамики приходится на второй и третий пласты (анимационный ролик 3) – здесь она встречена в некотором количестве в 27–29 разрезах из 59. Начиная с четвертого пласта, количество неопределенной керамики уменьшается, но она продолжает встречаться и глубже, вплоть до восьмого пласта, где она найдена в единственном разрезе В-6 (21 гр.). Максимальное количество неопределенной керамики – чуть более 60 гр. – найдено в третьем пласте разреза В-10. В большинстве разрезов на один пласт приходится не более 20 гр. подобной керамики.

4. Заключение. Полученные результаты геостатистического анализа распределения керамики на исследуемой почвенными разрезами территории в окрестностях укрепления Подкумское 3 позволили сделать следующие наблюдения. Проведенный анализ показал разные зоны сельскохозяйственного освоения изучаемой территории в кобанскую и аланскую эпоху. Так, для эпохи позднего бронзового – раннего железного века характерна максимальная концентрация находок керамики в зоне разреза А-8, где, вероятно, находилось поселение кобанского времени. Ресурсная зона данного поселения, на которую с видимой интенсивностью вносились органические удобрения, охватывала пространство на удалении в 200–300 м от данного места обитания. Кроме того, значительное количество фрагментов посуды этого периода, а также костей животных были обнаружены в разрезе D-28, устроенном на земледельческой террасе, в юго-западной части изучаемого комплекса. Расположенные поблизости на других террасах разрезы также несли значительное количество керамического материала и маркировали зону сельскохозяйственного освоения уже другого поселения, которое, скорее всего, находилось неподалеку, но не попало в изучаемое разрезами пространство.

Максимальная концентрация керамики эпохи раннего средневековья приходится на разрез А-3; далее наблюдаются постепенное снижение количества керамических фрагментов по мере удаления от укрепления Подкумское 3 и постепенное исчезновение аланской керамики в разрезах, находящихся на расстоянии более 650–700 м от крепости. Пространство интенсивного внесения удобрений в основном охватывало радиус около 500 м от раннесредневекового поселения, далее следовала зона менее интенсивной сельхозобработки. На расстоянии около 900 м в западном направлении имелась еще одна зона регулярной распашки (разрезы J-16 – N-16), очевидно, также относящаяся к укреплению Подкумское 3. Далее после некоторого перерыва в юго-западной части исследуемой в 2018 г. территории обнаружена периферия ресурсной зоны следующего раннесредневекового поселения, которое также не попало в сферу нашего внимания.

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

Основной вывод, полученный при геостатистическом анализе керамики из почвенных разрезов, может сводиться к следующему наблюдению над ее количеством: более 100 гр. керамики на один пласт в разрезе, скорее всего, отражает присутствие поселенческих структур и культурных слоев, либо соотносится с погребенной почвой земледельческих террас кобанского периода, где насыщенность керамическим материалом достаточно высока [Борисов, Коробов, 2013. С. 184–185]. Зона интенсивного внесения удобрений, вероятно, совпадает с количеством керамики более 30 гр. на почвенный разрез в целом, тогда как менее интенсивная обработка маркируется меньшим количеством керамики. Данные показатели будут соотноситься впоследствии с результатами химического и микробиологического анализа почвенных свойств изученных разрезов.

Работы по изучению ресурсной зоны разновременных поселений Кисловодской котловины будут продолжены в следующем полевом сезоне.

 

Литература

Балановская Е.В., Балановский О.П., 2007. Русский генофонд на Русской равнине. М.: ООО «Луч». 416 с.  назад

Борисов А.В., Коробов Д.С., 2013. Древнее и средневековое земледелие в Кисловодской котловине: итоги почвенно-археологических исследований. М.: ТАУС. 272 с.  назад

Демьянов В., Савельева Е., 2010. Геостатистика. Теория и практика. М.: Наука. 327 с.  назад

Ковалевская В.Б., 2016. Кавказ в свете археологии и геногеографии // РА. № 3. С. 45–59.  назад

Коробов Д.С., 2017а. Новые результаты компьютерного ГИС-моделирования ресурсных зон аланских поселений Кисловодской котловины I тыс. н.э. [Электронный ресурс] // Археология и геоинформатика. Вып. 8 / Отв. ред.: Г.Е. Афанасьев, Д.С. Коробов. М.: ИА РАН. DVD-ROM.  назад

Коробов Д.С., 2017б. Система расселения алан Центрального Предкавказья в I тыс. н.э. (ландшафтная археология Кисловодской котловины). Т. 1. М.; СПб.: Нестор-История. 384 с.  назад

Коробов Д.С., 2017в. Система расселения алан Центрального Предкавказья в I тыс. н.э. Т. 2. Каталог поселений Кисловодской котловины. М.; СПб.: Нестор-История. 312 с.  назад

Мыслыва Т.Н., Куцаева О.А., Белявский Ю.А., 2017а. Эффективность методов интерполяции на основе ГИС для оценки пространственного распределения рН в почве // Вестник Харьковского национального аграрного университета. Серия: Почвоведение, агрохимия, земледелие, лесное хозяйство, экология почв. № 2. С. 293–302.  назад

Мыслыва Т.Н., Куцаева О.А., Подлесный А.А., 2017б. Сравнение эффективности методов интерполяции на основе ГИС для оценки пространственного распределения гумуса в почве // Вестник Белорусской государственной сельскохозяйственной академии. № 4. С. 146–152.  назад

Чернышева Е.В., Борисов А.В., Коробов Д.С., 2016. Биологическая память почв и культурных слоев археологических памятников. М.: ГЕОС. 240 с.  назад

Barker G.W.W., 1975. Prehistoric territories and economies in Central Italy // Paleoeconomy / Ed. by E.S. Higgs. Cambridge: Cambridge University Press. P. 111–175.  назад

Bintliff J., 2000. The concepts of ‘site’ and ‘offsite’ archaeology in surface artefact survey // Non-Destructive Techniques Applied to Landscape Archaeology / Ed. by M. Pasquinucci, F. Trément. Oxford: Oxbow Books. P. 200–215. (The Archaeology of Mediterranean Landscapes; vol. 4).  назад

Chernysheva E.V., Korobov D.S., Khomutova T.E., Borisov A.V., 2015. Urease activity in culture layers at archaeological sites // Journal of Archaeological Science. Vol. 57. P. 24–31.  назад

Chernysheva E.V., Korobov D.S., Borisov A.V., 2017. Thermophilic microorganisms in arable land around medieval archaeological sites in Northern Caucasus, Russia: Novel evidence of past manuring practices // Geoarchaeology. Vol. 32, iss. 4. P. 494–501.  назад

Conolly J., Lake M., 2006. Geographical Information Systems in Archaeology. Cambridge: Cambridge University Press. 358 p.  назад

O'Connor T., Evans J.G., 2005. Environmental Archaeology: Principles and Methods. Stroud: Sutton Publishing. 256 p.  назад

Foley R., 1977. Space and Energy: A Method for Analysing Habitat Value and Utilization in Relation to Archaeological Sites // Spatial Archaeology / Ed. by D. Clarke. London; New York; San Francisco: Academic Press. P. 163–187.  назад

Ford S., Bowden M., Gaffney V., Mees G.C., 1994. The “Celtic” Field Systems on the Berkshire Downs, England // The Archaeology of Garden and Field / Ed. by N.F. Miller, K.L. Gleason. Philadelphia: University of Pennsylvania Press. P. 153–167.  назад

Higgs E.S., Jarman M.R., 1975. Paleoeconomy // Paleoeconomy / Ed. by E.S. Higgs. Cambridge: Cambridge University Press. P. 1–8.  назад

Hodder I., Orton C., 1976. Spatial Analysis in Archaeology. Cambridge; New York: Cambridge University Press. 278 p. (New studies in archaeology; 1).  назад

Jarman M.R., 1972. A territorial model for archaeology: a behavioral and geographical approach // Models in Archaeology / Ed. by D.L. Clarke. London: Methuen & Co Ltd. P. 705–733.  назад

Jarman M.R., Vita-Finzi C., Higgs E.S., 1972. Site catchment analysis in archaeology // Man, settlement and urbanism / Ed. by P.J. Ucko, R. Tringham, G.W. Dimbleby. London: Duckworth. P. 61–66.  назад

Korobov D.S., Borisov A.V., 2013. The origins of terraced field agriculture in the Caucasus: new discoveries in the Kislovodsk basin // Antiquity. Vol. 87, iss. 338. P. 1086–1103.  назад

Lloyd C.D., Atkinson P.M., 2004. Archaeology and geostatistics // Journal of Archaeological Science. Vol. 31 , iss. 2. P. 151–165.  назад

Miller N.F., Gleason K.L., 1994. Fertilizer in the Identification and Analysis of Cultivated Soil // The Archaeology of Garden and Field / Ed. by N.F. Miller, K.L. Gleason. Philadelphia: University of Pennsylvania Press. P. 25–43.  назад

Sjöberg A., 1976. Phosphate Analysis of Anthropic Soils // Journal of Field Archaeology. Vol. 3, no. 4. P. 447–454.  назад

Thurston T.L., 2001. Landscapes of Power, Landscapes of Conflict. State Formation in the South Scandinavian Iron Age. New York; Boston; Dordrecht: Kluwer Academic Publishers. 344 p.  назад

Wheatley D., Gillings M., 2002. Spatial Technology and Archaeology. The archaeological applications of GIS. London; New York: Taylor & Francis. 283 p.  назад

Widgren M., 1983. Settlement and farming systems in the early Iron Age. A study of fossil agrarian landscapes in Östergötland, Sweden. Stockholm: Almquist & Wiksell Int. 132 p. (Acta Universitatis Stockholmensis, Stockholm Studies in Human Geography; 3).  назад

Wilkinson T.J., 1989. Extensive Sherd Scatters and Land-Use Intensity: Some Recent Results // Journal of Field Archaeology. Vol. 16, no. 1. P. 31–46.  назад

Williamson T.M., 1984. The Roman Countryside: Settlement and Agriculture in N. W. Essex // Britannia. Vol. 15. Р. 225–230.  назад