REMOTE SENSING

Earth Engine & Geemap: 종 분포 모델링(SDM) 구현 (4) 모델 적합 및 예측

유병혁 2023. 8. 11. 20:50

안녕하세요? 이번 글은 시리즈 글의 일부로 Google Earth Engine geemap을 이용한 종 분포 모델링(SDM: Species Distribution Modeling) 구현 방법을 소개해 보겠습니다. 이 글의 내용은 스미소니언 보전생물연구소 연구진 분들이 제공한 JavaScript 소스코드를 Python으로 변환하여 수정, 보완한 것입니다.

 

Crego, R. D., Stabach, J. A., & Connette, G. (2022). Implementation of species distribution models in Google Earth Engine. Diversity and Distributions, 28, 904916. https://doi.org/10.1111/ddi.13491

 

사례 연구 1: 출현-전용 데이터를 이용한 팔색조(Pitta nympha) 서식지 적합성 및 예측 분포 모델링

5. 모델 적합

이제 모델을 적합(fit)할 수 있습니다. 모델을 적합한다는 것은 데이터의 패턴을 파악하여 모델의 매개변수(가중치와 편향)를 조정하는 과정을 의미합니다. 이렇게 하면 모델은 새로운 데이터를 입력 받았을 때 더 정확한 예측을 수행할 수 있게 됩니다. 일단, 모델을 적합하는 SDM() 함수를 정의합니다.

 

모델 훈련 및 검증 데이터 분할을 위해 공간 블록 교차-검증 기법(spatial block cross-validation technique)이 구현되었습니다. 사용자는 무작위 블록 분할로 여러 반복을 실행하여 훈련 및 검증 데이터셋을 실행할 수 있습니다.

그림 출처: https://onlinelibrary.wiley.com/doi/10.1111/ddi.13491

GEE에서 구현 가능한 몇 가지 비모수적 분류 알고리즘(non-parametric classification algorithms)이 있습니다. 비모수적 분류 알고리즘은, 데이터의 분포나 형태에 대한 가정을 하지 않고 주어진 데이터로부터 패턴을 추출하여 분류 작업을 수행하는 알고리즘을 말합니다. 랜덤 포레스트(random forest), 서포트 벡터 머신(support vector machine), 분류 및 회귀 트리(classification and regression trees), 최대 엔트로피(maximum entropy), 그래디언트 부스팅(gradient boosting)이 포함됩니다. 여기서는 랜덤 포레스트(random forest)를 사용해 봅니다.

  • numberOfTrees: 생성할 결정 트리의 수
  • variablesPerSplit: 분할당 변수의 수. 기본값: 변수 수의 제곱근 사용
  • minLeafPopulation: 훈련 셋에 최소한 이 개수의 지점이 포함된 노드만 생성
  • bagFraction: 트리당 백에 대한 입력 비율
  • maxNodes: 각 트리의 최대 리프 노드 수. 기본값: 제한 없음
  • seed: 무작위화 시드

표 출처: https://developers.google.com/earth-engine/apidocs/ee-classifier-smilerandomforest

임의-비출현 데이터 수는 출현 데이터 수와 균형을 이루는데 이 방법이 머신러닝 분류기에서 권장되기 때문입니다. setOutputMode() 함수를 사용하여 PROBABILITY(출현 확률) 및 CLASSIFICATION(이진 출현/비출현 지도)으로 결과를 얻습니다.

def SDM(x):
    Seed = ee.Number(x)

    # 훈련 및 검증을 위한 무작위 블록 분할
    GRID = ee.FeatureCollection(Grid).randomColumn(seed=Seed).sort('random')
    TrainingGrid = GRID.filter(ee.Filter.lt('random', split)) # 훈련용 격자
    TestingGrid = GRID.filter(ee.Filter.gte('random', split)) # 시험용 격자
    
    # 출현 지점
    PresencePoints = ee.FeatureCollection(Data)
    PresencePoints = PresencePoints.map(lambda feature: feature.set('PresAbs', 1))
    TrPresencePoints = PresencePoints.filter(ee.Filter.bounds(TrainingGrid)) # 훈련용 출현 지점
    TePresencePoints = PresencePoints.filter(ee.Filter.bounds(TestingGrid)) # 검증용 출현 지점
    
    # 임의-비출현 지점
    TrPseudoAbsPoints = AreaForPA.sample(region=TrainingGrid,
                                         scale=GrainSize,
                                         numPixels=TrPresencePoints.size().add(300),
                                         seed=Seed,
                                         geometries=True)
    # 훈련용 출현 지점과 동일한 수의 임의-비출현 지점
    TrPseudoAbsPoints = TrPseudoAbsPoints.randomColumn().sort('random').limit(ee.Number(TrPresencePoints.size()))
    TrPseudoAbsPoints = TrPseudoAbsPoints.map(lambda feature: feature.set('PresAbs', 0))
    
    TePseudoAbsPoints = AreaForPA.sample(region=TestingGrid,
                                         scale=GrainSize,
                                         numPixels=TePresencePoints.size().add(100),
                                         seed=Seed,
                                         geometries=True)
    # 검증용 출현 지점과 동일한 수의 임의-비출현 지점
    TePseudoAbsPoints = TePseudoAbsPoints.randomColumn().sort('random').limit(ee.Number(TePresencePoints.size()))
    TePseudoAbsPoints = TePseudoAbsPoints.map(lambda feature: feature.set('PresAbs', 0))

    # 훈련 및 임의-비출현 지점 병합
    trainingPartition = TrPresencePoints.merge(TrPseudoAbsPoints) # 훈련용 지점
    testingPartition = TePresencePoints.merge(TePseudoAbsPoints) # 시험용 지점

    # 훈련용 지점에서 예측 변수 이미지의 공변량 값 추출 
    trainPixelVals = predictors.sampleRegions(collection=trainingPartition,
                                              properties=['PresAbs'],
                                              scale=GrainSize,
                                              tileScale=16,
                                              geometries=True)

    # 랜덤 포레스트 분류기(Random Forest classifier)
    Classifier = ee.Classifier.smileRandomForest(
        numberOfTrees=500, # 생성할 결정 트리의 수
        variablesPerSplit=None, # 분할당 변수의 수. 기본값: 변수 수의 제곱근 사용
        minLeafPopulation=10, # 훈련 셋에 최소한 이 개수의 지점이 포함된 노드만 생성
        bagFraction=0.5, # 트리당 백에 대한 입력 비율
        maxNodes=None, # 각 트리의 최대 리프 노드 수. 기본값: 제한 없음
        seed=Seed # 무작위화 시드
    )
    
    # 출현 확률(Presence probability)
    ClassifierPr = Classifier.setOutputMode('PROBABILITY').train(trainPixelVals, 'PresAbs', bands)
    ClassifiedImgPr = predictors.select(bands).classify(ClassifierPr)
    
    # 이진 출현/비출현 지도(Binary presence/absence map)
    ClassifierBin = Classifier.setOutputMode('CLASSIFICATION').train(trainPixelVals, 'PresAbs', bands)
    ClassifiedImgBin = predictors.select(bands).classify(ClassifierBin)
   
    return ee.List([ClassifiedImgPr, ClassifiedImgBin, trainingPartition, testingPartition])

공간 블록은 모델 적합의 경우 70%, 모델 검증의 경우 30%로 각각 분할됩니다. 임의-비출현 데이터는 반복할 때마다 각 훈련 및 검증 셋 내에서 무작위로 생성됩니다. 결과적으로 매 실행 때마다 모델 적합 및 검증을 위한 서로 다른 출현 및 임의-비출현 셋을 갖게 됩니다.

split = 0.70  # 공간 블록은 모델 적합의 경우 70%, 모델 검증의 경우 30%로 각각 분할
numiter = 10 # 반복 횟수

runif = lambda length: [random.randint(1, 1000) for _ in range(length)] # 시드
results = ee.List(runif(numiter)).map(SDM).flatten()

5. 모델 예측

이제 모델 예측을 추출하여 지도에 표시할 수 있습니다. 서식지 적합성 지도(Habitat suitability map)잠재 분포 지도(Potential distribution map)를 가시화한 결과는 다음과 같습니다.

# 서식지 적합성 지도(Habitat suitability map)
images = ee.List.sequence(0, ee.Number(numiter).multiply(4).subtract(1), 4).map(lambda x: results.get(x))
ModelAverage = ee.ImageCollection.fromImages(images).mean()

Map = geemap.Map()
vis_params = {
    'min': 0,
    'max': 1,
    'palette': ['#440154', '#482677', '#404788', '#33638D', '#287D8E',
          '#1F968B', '#29AF7F', '#55C667', '#95D840', '#DCE319']}
Map = geemap.Map()
Map.addLayer(ModelAverage, vis_params, 'Habitat suitability')
Map.add_colorbar(vis_params, label="Habitat suitability", orientation="horizontal", layer_name="Habitat suitability")
Map.addLayer(Data, {'color':'red'}, 'Presence')
Map.centerObject(Data.geometry(), 7)
Map

# 잠재 분포 지도(Potential distribution map)
images2 = ee.List.sequence(1, ee.Number(numiter).multiply(4).subtract(1), 4).map(lambda x: results.get(x))
DistributionMap = ee.ImageCollection.fromImages(images2).mode()

vis_params = {
    'min': 0,
    'max': 1,
    'palette': ['white', 'green']}
Map = geemap.Map()
Map.addLayer(DistributionMap, vis_params, 'Potential distribution')
Map.addLayer(Data, {'color':'red'}, 'Presence')
Map.add_colorbar(vis_params, label="Potential distribution", discrete=True, orientation="horizontal", layer_name="Potential distribution")
Map.centerObject(Data.geometry(), 7)
Map

1920_ee_geemap_종 분포 모델링(SDM) 구현 (4) 모델 적합 및 예측.ipynb
0.89MB