Mimika Land Rights Project – West Papua
Prepared by Mike Alonzo (GISCorps Volunteer) and Nabil Ahmed (Forensic Architecture)
West Papua is the remote Melanesian province of Indonesia that has an abundance of gas, mineral and timber resources (Fig 1). Together with Papua New Guinea, in terms of biodiversity, West Papua contains 5% of the world’s tropical high-biodiversity terrestrial ecosystem. At the same time West Papua is a conflict zone with heavy Indonesian military presence, as there is a decade’s long independence movement waged by Papuans. The military has a long history of dealing heavy-handed with political protest and civilians alike and there are numerous cases of well-documented human rights violations that amounts to a slow genocide.
Figure 1: The study area in Western Papua
West Papua remains a region under close supervision by the military and a difficult area to operate in for human rights workers, environmental activists and journalists. The Grasberg mine operated by Freeport-Rio Tinto is located in Amungme territory in the Carstenz range. It is considered to hold the planet’s largest combined reserve of copper and gold. Copper from the Grasberg mine is the cheapest in the world.
The Grasberg mine has been dumping tailings into the Ajkwa River deposition area downstream in Komoro territory since operation in the 1980’s and is the site of an environmental disaster. At the same time, Freeport-Rio Tinto have a documented history of financing Indonesian military presence in West Papua, which activists claim is supporting a slow genocide.
Nabil Ahmed (Forensic Architecture and Earth Sensing Association) has launched an investigation into documenting the Grasberg mine environmental disaster on behalf of Octavianus Mote, Papuan negotiator and in collaboration with TAPOL, the human rights campaign organization for Indonesia based in the UK. The spatial analysis conducted by GISCorps Remote Sensing Specialist is a crucial preliminary element of the investigation, especially for an area that is extremely difficult to reach as a zone of conflict.
Objectives and Deliverables (i.e., tasks for the GISCorps Remote Sensing Specialist)
Objective #1: Map vegetation disturbance in Western Papua that can be attributed to either acid mine drainage from the Freeport-Rio Tinto mine or to population pressures (clearing for agriculture, urbanization) associated with increased mining activity.
- Deliverable 1: Animated visualization of land cover change from 1987 – 2014.
- Deliverable 2: Dated disturbance map from 1987 – 2014.
- Deliverable 3: Landscape complexity map illustrating the number of “states” of every pixel from 1987 – 2014.
Objective #2: Quantify the growth of the acid mine drainage area (ADA) by year and correlate with existing mine production data.
- Deliverable 1: Amount of land disturbed by year
- Deliverable 2: Visual comparison of disturbance with mine tailing production data
- Western Papua is a challenging area for remote sensing analysis because it is extremely cloudy. Also, there is very little access to high resolution imagery or in situ data for validating remote sensing analyses.
- Given that our goal was to evaluate change over time at moderate resolution with a focus on the 1990’s (the decade where the Freeport-Rio Tinto mining operation was in heaviest production), we opted for Landsat imagery. The additional advantages of Landsat include: (1) ease of access to free reflectance (atmospherically corrected) imagery) and (2) 16-day revisit cycle gave us a better chance of finding cloud-free pixels at some point for every year.
- Specifically we downloaded 199 LEDAPS surface reflectance images from Landsat 4/5 TM and Landsat 7 ETM+. This was basically every image that was remotely cloud free (~20% cloud free) over the full time period (1987 – 2014).
- Ancillary data was mostly from a report titled: “The Environmental Impacts of Freeport-Rio Tinto’s Copper and Gold Mining Operation in
- Papua ” by WALHI (Friends of the Earth Indonesia)
- All data processing was completed in Matlab. Code will be made available by Mike Alonzo after publication of an associated journal article.
- Overview: The basic plan was to create an algorithm which could track land cover change (anywhere in the world, potentially) based on change in any spectral index including but not limited to the Normalized Difference Vegetation Index (NDVI). The work described in this report was heavily inspired by the USDA Forest Service’s Landtrendr project. An adapted algorithm was required, however, because of the extreme cloudiness of the Western Papua region.
- All cloudy areas were masked using the cfmask band (which only masks the most obvious clouds) and the 199 images were clipped to our study area extent. Significant clouds and haze remained (Fig 2).
Figure 2: Study area during a particularly cloudy year. Solid green areas are the cfmask but note how many additional clouds remain unmasked.
- All reflectance images were transformed to NDVI which responds to the presence and density of vegetation. NDVI is calculated as (Near Infrared – Red) / (Near infrared + Red). We now have a data cube that is 2000 pixels (60 km north – south), 1300 pixels (39 km east – west) and 199 pixels in depth (1 for each date) (Fig 3).
- If one clicks on a pixel that falls in the ADA, wanting to see the NDVI values over time, it will look something like this (Fig 4):
- If one clicks on stable forest, the NDVI values will look something like this (Fig 5):
Figure 5: NDVI change for growing, undisturbed forest. The many instances of low NDVI can be considered noise related to atmospheric contamination (mostly clouds).
- The goal now is to develop an algorithm that can fit several lines in a piece-wise fashion to this rather noisy data in order to show: (1) the number of “states” of this pixel (e.g., 3 states could be stable forest -> ADA disturbance event -> stable ADA); (2) the type of disturbance (e.g., urbanization vs. ADA); and (3) the date of the disturbance. All of this must happen while avoiding “distraction” from the many data points that are bad because of cloud/haze intrusion.
- The basis for the chosen algorithm is Matlab code called Shape Language Modeling (SLM) which was written by John D’Errico. This code makes use of non-linear optimization to simultaneously estimate the dates of a state change (e.g., change from stability to decline) as well as the associated NDVI value at that point. Figure 4 shows the result of this code visually using the red line for ADA disturbance. Figure 5 shows the single line for undisturbed forest. Figure 6 shows the process of urbanization. This is different from the ADA disturbance because it is more gradual (i.e., the slope of the red line is not as steep).
Figure 6: Urbanization characterized by a somewhat gradual loss of vegetation followed by a period of stability.
- We ran this code for all 2.6 million pixels in the study area. This is very processing intensive and could literally take a month to run on a single processor. Luckily, many folks in the remote sensing community have access to some kind of multicore server. For this analysis, we used 32 processors in parallel and completed the study area in approximately 1.5 days.
- Disturbance dates were assigned based on a change from higher than 0.40 NDVI to below 0.40 NDVI and a slope of the cross-over segment that was adequately high (a means for separating the somewhat ill-defined processes of “disturbance” and “decline”).
- The year-by-year increase in the ADA was simply calculated as the number of pixels disturbed per year. This does include urbanized pixels as well but that number is fairly small and could be separated out simply if required. One pixel = 900 m2.
- The land cover change animation: This deliverable was dealt with separately. Here we wanted to produce one relatively cloud-free composite image per year to add to an animation. To this end, we selected the pixel-date in a given year for each pixel location. “Best” was defined as smallest band1:band4 (blue to near infrared) ratio because high blue reflectance is broadly indicative of atmospheric contamination.
- Land cover change animation (Click here to download the animation):
- Some of the typical “pixel trajectories” (aka land cover histories) are shown in Figures 4, 5, and 6 above.
- Here is the map of estimated disturbance date. Black areas were not disturbed (Fig 7).
- Here is the map of landscape complexity. The interpretations of the complexities (with varying degrees of confidence are):
1 = stable forest (or less likely, stable river)
2 = decline/disturbance à stable ADA or urban area
3 = stable forest à disturbance à stable ADA4 and 5: Probably complex, braided river dynamics possibly including disturbance and vegetation regrowth. (Fig 8)
Figure 8: Landscape complexity. Or, how many things have happened to a given pixel. Most of the ADA has complexity of 3 which means it was first forested, then disturbed, and is now “river”.
- Amount of land disturbed by year overlaid on the WALHI figure for mine daily mine tailings production. (Fig 9)
Figure 9: Background (and left axis) is data from WALHI report on tailings production per year. Foreground (and right axis) is yearly growth of ADA. We hypothesize that after about 2000, the tailings were simply reaching the sea, thus the slowed growth of the ADA. This has not been rigorously confirmed.
- Our findings clearly show quantifiable environmental damage to the Ajkwa River ADA that can be attributed to mining activities over the full time period (1987 – 2014)
- It further makes clear that as the mining activities increased over the years so did the tailing production and environmental disturbance as a consequence.
- The report indicates further research is required to document more fine-grained analysis in relation to determining vegetation loss by tree and forest type and water quality that is the source of livelihood and subsistence by local Komoro communities. It is important for us to connect the remote sensing work to analyze and evidence socio-environmental violations and conflict.
- The findings of the report will help in designing future ground truth data collection.
- The analysis and results of this report rely on what is likely the most sophisticated remote sensing study of this area. We will make available any outputs that may be of use to similar academic (focused on the same area or otherwise) or GISCorps projects.