Способ актуализации лесотаксационных данных на основе машинного обучения

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

 

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

Известен (RU, патент 2080051, опубл. 09.05.1995) способ определения запаса насаждений. Согласно известному способу получают изображение лесного массива, разбивают изображение на мозаику участков, последовательно преобразуют изображение каждого участка в матрицу цифровых отсчетов пространственной зависимости яркости I (х, у) размерностью I m х m I элементов, вычисляют характеристики электрического сигнала матрицы: огибающую пространственного спектра, среднеквадратическое отклонение, рассчитывают запас насаждения по статистическим зависимостям параметров сигнала и таксационных характеристик.

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

Известен также (RU, патент 2133565, опубл. 27.07.1999) способ оценки запаса насаждений путем получения изображения лесного массива, преобразования его в дискретную матрицу отсчетов яркости, нахождения всех локальных максимумов матрицы, расчета среднего расстояния между максимумами, выделения контуров вокруг каждого максимума методами пространственного дифференцирования и оценки запаса по аналитическим зависимостям между рассчитанными числовыми характеристиками матрицы и элементами дерева среднего класса Лорея.

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

Также известен (RU, патент 2038001, опубл. 27.06.1995) Способ оценки состояния лесов. В известном способе осуществляют дистанционные спектрометрические измерения, с борта орбитальной станции, коэффициентов спектральной яркости зондируемой лесной площади в синем В, зеленом G и красном R участках видимого спектра, рассчитывают значения хроматических коэффициентов жизненности

g=G/(B+Gh-R)

и красного поражения

r=R/(B+G+R),

вычисляют регрессию хроматических коэффициентов

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

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

Известен также (RU, патент 2406295, опубл. 20.12.2010) способ экологического мониторинга лесов, включающий дистанционную регистрацию полей яркости лесной растительности аэрокосмическими средствами, причем дистанционную регистрацию полей яркости лесной растительности осуществляют путем зондирования много- или гиперспектральным датчиком в зеленой G (450-550 нм), красной R (550-670 нм) и ближней инфракрасной БИК (670-950 нм) зонах спектра с одновременным получением цифровых изображений для каждой зоны и проводят математическую обработку полученной информации (патент RU №2406295). Оперативность, производительность и глобальность контроля состояния лесных массивов обеспечивается особенностями космической съемки.

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

Известный источник информации принят в качестве ближайшего аналога.

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

Технический результат, достигаемый при реализации разработанного способа, состоит ускорении работ по инвентаризации лесов.

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

В некоторых вариантах реализации разработанного способа для проведения таксации территории используются данные дистанционного зондирования с разрешением 10-30 м по не менее, чем 4 спектральным каналам в диапазоне 0,4-13 мкм за временной период, охватывающий не менее, 5 периодов наблюдения, расположенными в различных фазах вегетационного периода, охватывающие раннюю весну, позднюю весну, середину лета, конец лета и осенний период.

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

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

Для решения данной технической задачи необходимо произвести подготовку данных дистанционного зондирования и подготовку существующих инвентаризационных данных.

Используемая дистанционная информация делится на две основных группы:

1) мультиспектральная спутниковая информация и ее производные,

2) цифровая модель рельефа и ее производные.

Мультиспектральная спутниковая информация есть результат измерения структуры потока отраженной солнечной радиации. Отражение солнечной радиации можно рассматривать как некоторый эмерджентный эффект функционирования каждой элементарной единицы земной поверхности. Тем самым формально мультиспектральное отражение содержит информацию о априори неизвестных процессах трансформации энергии ландшафтным покровом. Наиболее подходящими материалами дистанционного зондирования таксируемой территории являются системы наблюдения с разрешением 10-30 м по четырем-десяти спектральным каналам в диапазоне 0,4-13 мкм. Мультиспектральное отражение существенно меняется по сезонам года, причем эти изменения различны у различных экосистем. Экосистемы не различимые по спектру в один сезон года могут хорошо различаться в другой. Так, например, хвойные леса от лиственных наилучшим образом можно различить по NDVI (Normalized Difference Vegetationlndex - нормализованный относительный индекс вегетации) зимой. Это определяет необходимость использовать для построений карт ландшафтного покрова и отдельных его компонентов спектральные съемки за различные сроки наблюдения. Необходимо составить базу данных дистанционной информации, описывающую ход функционирования вегетации в течении вегетационного сезона с не менее чем 5 периодами наблюдения расположенными в различных месяцах, охватывающие (в контексте климатических условий) раннюю весну, позднюю весну, середину лета, конец лета и осенний период.

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

Расчет термодинамических характеристик поверхности отраженных в мультиспектральной мультиспектральной съемке для единицы поверхности осуществляется с использованием оценки дистанции между распределением мощностей по спектру поглощенной солнечной энергии и равновесным состоянием - гипотетическим поглощением солнечной энергии пропорционально распределению мощностей в спектре солнечной постоянной. Степень отклонения спектра отраженной энергии от спектра приходящей измеряют энтропией (информацией) Кульбака (Assessing the thermodynamic variables of landscapes in the southwest part of east european plain in russia using the modis multispectral band measurements / Y.G. Puzachenko, R.B. Sandlersky, A.N. Krenke, A. Olchev // Ecological Modelling. - 2016. - Vol.319. - P. 255-274). Так как поглощающей поверхностью является в основном растительный покров, то изменения морфологической и биохимической структуры во времени и в пространстве могут возникать только в результате извлечения информации из среды, в том числе и из структуры солнечного спектра. Приращение информации по Кульбаку равно нулю, если распределение поступающей и отраженной радиации по полосам спектра тождественно и, соответственно приемник информации равновесен передатчику. Если информация Кульбака больше нуля, то можно говорить о приращении информации в приемнике и отражающая поверхность неравновесна относительно спектра солнечной радиации.

Приращение информации по Кульбаку (К, нит) рассчитывается по формуле (там же):

где - доля пришедшей энергии в спектральном диапазоне ν от суммарной поступившей энергии; - доля отраженной энергии в спектральном диапазоне ν от суммарной отраженной энергии, - пришедшая энергия в спектральном диапазоне ν, - отраженная энергия в спектральном диапазоне ν.

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

Полезную работу солнечной радиации (Ех) рассчитывают как:

где Eⁱⁿ - приходящая солнечная энергия, E^out - отраженная солнечная энергия, R=Eⁱⁿ-E^out - поглощенная энергия, a E^out/Eⁱⁿ - альбедо.

При фиксированном альбедо, чем больше приращение информации, тем больше эксэргия солнечной радиации. Обратим внимание на то, что вклад первого члена уравнения в эксэргию по абсолютной величине относительно не велик.

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

Далее рассчитывают связанная энергия (STW):

STW=TW*S_out,

где TW - тепловой поток, от деятельной поверхности, зафиксированный 6-м, длинноволновым каналом.

Приращение внутренней энергии системы (DU) - переход поглощенной солнечной энергии во внутреннюю энергию системы, оценивают как остаток в балансовом уравнении поглощенной энергии (R):

DU=R-Ex-STW.

Для оценки непосредственных затрат энергии на производство биологической продукции используют разность отраженной энергии в красном (RED) и ближнем инфракрасном диапазонах (NIR) - вегетационный индекс (VI), или, иначе индекс продуктивности. Расчет VI базируется на двух наиболее стабильных (не зависящих от прочих факторов) участках спектральной кривой отражения сосудистых растений. В красной области спектра (0.60-0.70 мкм) лежит максимум поглощения солнечной радиации хлорофиллом высших сосудистых растений, а в инфракрасной области (0.70-1.0 мкм) находится область максимального отражения клеточных структур листа. Высокая фотосинтетическая активность ведет к меньшему отражению в красной области спектра и большему в ближней инфракрасной. Здоровая плотная растительность поглощает большую часть излучения в красной области спектра и отражает большую часть в инфракрасной, соответственно, с уменьшением плотности, доля отраженного излучения в красной области увеличивается, а в инфракрасной уменьшается. Индекс рассчитывается как:

где - отражение солнечной энергии в диапазоне 0.77-0.90 мкм, - отражение солнечной энергии в диапазоне 0.63-0.69 мкм.

В итоге рассчитывают следующие термодинамические характеристики:

Eⁱⁿ - приходящая солнечная энергии, Вт/м²;

E^out - отраженная солнечная энергия Вт/м²;

R - поглощенная энергия, Вт/м²;

К - приращение информации (энтропия Кульбака), нит;

S_out - энтропия отраженной солнечной радиации, нит;

TW - тепловой поток от деятельной поверхности, Вт/м²;

Т - температура деятельной поверхности, °С;

Ех - эксергия солнечной радиации, Вт/м²;

STW - связанная энергия, Вт/м²;

DU - приращение внутренней энергии, Вт/м²;

VI - индекс продуктивности, Вт/м².

Цифровую модель рельефа подвергают процедуре иерархического анализа ведущих компонент его организации на основе двухмерного спектрального анализа. Последовательно решают следующие частные задачи:

1. Расчет двухмерного спектра рельефа методом преобразования Хартли

2. Свертка двухмерного спектра в линейную форму по левому и правому секторам сдвига при оценки автокорреляционной функции

3. Оценка параметров регрессии «логарифм спектра - логарифм частоты» и определение на этой основе фрактальной размерности (D):

Log(Спектральная плотность)=а+b log(частота),

4. Вычисление фрактальной размерности: D=(7-b)/2

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

6. Вычисление остатков от линии регрессии и выделение иерархических структур внутри спектра методом итеративного моделирования остатков квадратичной функцией с выделением ее максимума.

7. Анализ спектра остатков и оценка возможного числа факторов, порождающих пространственную структуру рельефа.

8. Построение общей статистической модели иерархической организации рельефа.

Исходная лесотаксационная база данных представляет собой единообразный набор таблиц для каждого лесничества. Эти таблицы содержат характеристики индивидуальных выделов. Каждому выделу соответствует отдельная строка с его характеристиками и, уникальным для всей базы данных, номером. С использованием программных пакетов можно последовательно свести оцифровку и дистанционную информацию в одну ГИС. Исходные данные сводят в одну проекцию (как правило, соответствующую проекции дистанционных данных). Далее оцифровку лесотаксационного планшета преобразовывают в файл с расширением shape, а потом в grid, при этом для каждой ячейки сетки сохраняется индивидуальный идентификатор выдела, занимающего на ней наибольшую площадь. Размер ячеек растра задается в соответствии с принятым масштабом. При привязке растра часто необходимо выполнить его геометрическую коррекцию относительно снимка.

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

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

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

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

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

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

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

2. Способ по п. 1, отличающийся тем, что для проведения таксации территории используются данные дистанционного зондирования с разрешением 10-30 м по не менее чем четырем спектральным каналам в диапазоне 0,4-13 мкм за временной период, охватывающий не менее пяти периодов наблюдения, расположенных в различных фазах вегетационного периода, охватывающих раннюю весну, позднюю весну, середину лета, конец лета и осенний период.

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