Do parks that are more strictly protected have less deforestation?
In this problem, we examine the question of land protection status and its impacts on deforestation activity in the country of Malawi. We use the Hansen Global Forest Change dataset to determine changes in tree cover between 2000 and 2019.
Our analysis found that in 2019, 0.46% of deforestation occurred in less protected parks, while only 0.01% occurred within strictly protected areas. Our results therefore indicate that less protected parks in Malawi have more deforestation (Map 1). Still, it is important to consider the ways in which these results are shaped by the assumptions of our study. For instance, we based our categories of “strictly protected” and “less protected” areas on protection status level ascribed by the IUCN, and the amount of human activity allowed within such areas. However, we also included areas listed as “Not Reported” or “Not Applicable,” as well as “Pending,” in our “less protected” category. In Malawi, the majority of protected areas (90.98%) are of “Not Reported” Status, meaning we are drawing our conclusions from considerably vague information with regard to how areas are managed.
Another aspect that our analysis overlooks is the fact that ‘forest loss’ in 2019 counts as all pixels categorized as loss in the Hansen et al. (2013) dataset, rather than those categorized as loss in areas that exhibited tree cover >=30% in 2000, which was our threshold of 30% canopy cover needed to represent forested land. Hansen et al. (2013) define forest loss as “complete removal of tree cover canopy at the Landsat pixel scale [30 m2].” This means that some of the ‘loss’ in our study occurred in pixels that exhibited total tree cover of less than 30%, which was outside the realm of our study. This gave us an overestimation of loss in 2019. In order to correct this, we would need to add in a command that excludes loss pixels from 2019 that do not coincide with pixels with canopy cover >=30% from 2000.
Here is a code segment that can be run in Google Earth Engine to replicate our methods and results:
/* 1. MAKE VARIABLES FOR PARTS OF IMAGES WE WANT TO USE */
var treeCover = ee.Image("UMD/hansen/global_forest_change_2019_v1_7").select(['treecover2000']).gte(30); //call on the treecover band and select for pixels w value >=30% cover
var lossin2019 = ee.Image("UMD/hansen/global_forest_change_2019_v1_7").select(['lossyear']).eq(19); //call on the lossyear band; select for pixels of the year 2019
/* 2. SELECT PARKS IN MALAWI AND PLACE THEM INTO TWO GROUPS BASED ON STATUS */
var countries = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017'); //Call on the Country Boundaries dataset
var malawi = ee.Feature( //Create a variable that specifies the geometry of Malawi
countries
.filter(ee.Filter.eq('country_na', 'Malawi')) //Filter for the criteria corresponding to Malawi
.first()
);
var StrictPA = ee.FeatureCollection('WCMC/WDPA/current/polygons') //Call on the Protected Areas dataset
.filter(ee.Filter.and( //Create a system that uses an 'and' function (requires all listed criteria to be met)
ee.Filter.bounds(malawi.geometry()), //Filter by the boundaries of Malawi geometry defined above
ee.Filter.neq('IUCN_CAT', 'VI'), //Filter out categories of less protected areas using "notequal" function
ee.Filter.neq('IUCN_CAT', 'V'),
ee.Filter.neq('IUCN_CAT', 'Not Reported'),
ee.Filter.neq('IUCN_CAT', 'Not Applicable'),
ee.Filter.neq('STATUS', 'proposed')
))
.map(function(feature){
return malawi.intersection(feature); //Map the Intersect function to clip the intersecting areas by Malawi geometry to return only the protected area pixels that fall entirely within Malawi
});
var LessPA = ee.FeatureCollection('WCMC/WDPA/current/polygons') //Same as above ^ only .neq to criteria for Strictly Protected Areas
.filter(ee.Filter.and(
ee.Filter.bounds(malawi.geometry()),
ee.Filter.neq('IUCN_CAT', '1a'),
ee.Filter.neq('IUCN_CAT', '1b'),
ee.Filter.neq('IUCN_CAT', 'II'),
ee.Filter.neq('IUCN_CAT', 'III'),
ee.Filter.neq('IUCN_CAT', 'IV')
))
.map(function(feature){
return malawi.intersection(feature); //Map the Intersect function to clip the intersecting areas by Malawi geometry to return only the protected area pixels that fall entirely within Malawi
});
/* 2. CALCULATE TREE COVER IN STRICT AND LESS PA 2019 */
var statsStrict = lossin2019.reduceRegion({ //reduce the image by the geometry of Malawi --> creates a dictionary object (cannot be used in math functions)
reducer: ee.Reducer.sum(), //sum all pixels that fall within the specified geometry (px with value of 1 = loss for this band; else, 0)
geometry: StrictPA.geometry(), //set the geometry to the strictly protected areas
scale: 30, //px scale relative to real-world units
maxPixels: 1e9
});
/*
print('Loss in Strictly Protected Areas 2019:', //print the count information to the console (just to check)
statsStrict.get('lossyear') //specifies the band which contains the information
);
*/
var statsLess = lossin2019.reduceRegion({ //reduce the image by the geometry of Malawi
reducer: ee.Reducer.sum(), //sum all pixels w/in geometry (px with value of 1 = loss for this band; else, 0)
geometry: LessPA.geometry(), //specify less protected area geometry
scale: 30,
maxPixels: 1e9
});
/*
print('Loss in Less Protected Areas 2019:', //print the count information to the console (just to check)
statsLess.get('lossyear'));
*/
var tc2000malawi = treeCover.reduceRegion({ //reduce the region of the treecover2000 infor to just inside the bounds of Malawi geometry
reducer: ee.Reducer.sum(), //sum px inside (reducer is sum)
geometry: malawi.geometry(), //specify geometry of boundary
scale: 30,
maxPixels: 1e9
}); //returns a dictionary object--> need to convert to numeric before performing calculations
/*
print('Tree Cover in Less Protected Areas of Malawi (2000)', //print the count information to the console
tc2000malawi.get('treecover2000')); //call on particular band to be printed
*/
/* 4. CALCULATE PERCENTAGE LOSS IN 2019 */
print('Percent of Treecover Lost in Strictly Protected Areas (2019)', //name printed on console
ee.Number(statsStrict.get('lossyear')) //statsStrict --> numeric; specify lossyear band for calculations (specifies all px from strictly protected areas)
.divide(ee.Number(tc2000malawi.get('treecover2000'))) //divide statsStrict pixels by total Malawi treecover from 2000
.multiply(100), //multiply by 100 to make into a percent
('percent'));
print('Percent of Treecover Lost in Less Protected Areas (2019)',
ee.Number(statsLess.get('lossyear'))
.divide(ee.Number(tc2000malawi.get('treecover2000')))
.multiply(100), //^same as above, but calculations based on px in less protected areas
('percent'));
/* 5. VISUALIZE RESULTS */
Map.setCenter(34.3015, -13.2543, 8);
//center the map on Malawi and zoom
Map.setOptions('TERRAIN'); //Terrain as basemap
Map.addLayer(treeCover.clip(malawi).selfMask(), //treeCover layer from 2000, clipped to extent of Malawi, with a self mask to eliminate 0 value px from the picture
{palette: ['487552']}, //tree color
"Tree Cover in 2000"); //layer name
Map.addLayer(lossin2019.clip(malawi).selfMask(), //^same as above
{palette:'red'},
"Tree Cover Loss in 2019");
var blankImage = ee.Image().byte(); //creates a blank image onto which you can paint the lines
var LessPAOutlines = blankImage.paint({ //new container to store the boudaries of the countries in an image; .paint will put the lines in the image
featureCollection: LessPA, //name of vector data that you want to paint into the image
width: 1.5 //the size of the line that you want to paint (1 = thin)
});
Map.addLayer(LessPAOutlines, {palette: '38d1c6'}, "Less Protected Area Boundaries"); //Add boundary lines to map
var blankImage = ee.Image().byte(); //creates a blank image onto which you can paint the lines
var StrictPAOutlines = blankImage.paint({ //new container to store the boudaries of the countries in an image; .paint will put the lines in the image
featureCollection: StrictPA, //name of vector data that you want to paint into the image
width: 1.5 //the size of the line that you want to paint (1 = thin)
});
Map.addLayer(StrictPAOutlines, {palette: 'black'}, "Strictly Protected Area Boundaries"); //Add boundary lines to map
And a link to the code:
https://code.earthengine.google.com/e9d52892a61b89059fd6437c411bcd58
| Works Cited: M.C. Hansen, P.V. Potapov, R. Moore, M. Hancher, S.A. Turubanova, A. Tyukavina, D. Thau, S.V. Stehman, S.J. Goetz, T.R. Loveland, A. Kommareddy, A. Egorov, L. P. Chini, C.O. Justice, J. R.G. Townshend (2013). High-Resolution Global Maps of 21st-Century Forest Cover Change. Science, 342(6160), pp. 850-853. |