Estimating CO₂ Sequestration by Forests in Oita Prefecture, Japan, by Combining LANDSAT ETM+ and ALOS Satellite Remote Sensing Data

Abstract

CO₂ sequestration of the forests in Oita Prefecture, Japan, was estimated using satellite remote sensing data. First, hybrid classification of the optical LANDSAT ETM+ data was performed using GIS to produce a detailed land cover map. CO₂ sequestration for each forest type was calculated using the sequestration rates per unit area multiplied by the forest areas obtained from the land cover map This results in 3.57 MtCO₂/yr for coniferous, 0.77 MtCO₂/yr for deciduous broadleaf, and 2.25 MtCO₂/yr for evergreen broadleaf, equivalent to a total of 6.60 MtCO₂/yr for all the forest covers in Oita. Then, two different methodologies were used to improve these estimates by considering tree ages: the Normalized Difference Vegetation Index (NDVI) and the stem volume methods. Calculation using the NDVI method shows the limitation of this method in providing detailed estimations for trees older than 15 years, because of NDVI saturation beyond this age. In the stem volume method, tree ages were deduced from stem volume values obtained by using PALSAR backscattering data. Sequestration based on tree age forest subclasses yields 2.96 MtCO₂/yr (coniferous) and 0.31 MtCO₂/yr (deciduous broadleaf forests). These results show the importance of using not only detailed forest types, but also detailed tree age information for more realistic CO₂ sequestration estimates. In so doing, overestimation of the sequestration capacity of forests could be avoided, and the information on the status and location of forest resources could be improved, thereby leading to sounder decision making in sustainable management of forest resources.

1. Introduction

Global warming is known as one of the biggest environmental issues of the 21st century. Many studies have been made and conferences held, emphasizing the importance of this issue and the associated problems. World nations are working together on policies, not only for preventing the problem, but mostly for mitigating its impacts [1–3]. Japan is one of the top nations involved in activities related to this issue, and the Kyoto Protocol (adopted at the Kyoto Conference in 1997), provided new approaches to finding solutions [4], among which are the Clean Development Mechanism (CDM), Emission Trading, Joint Implementation (JI) and Carbon Sink. Carbon sink is the focal point in the present research context, as 68.9% of land in Japan are covered by forests [5]. Japan is considering the forests to be the top candidate for its CO₂ reduction process [6]. Many prefectures in Japan have been researching carbon sequestration by forests and have come up with various policies for reducing CO₂. However, some prefectures have not yet been involved in such processes focusing on CO₂ reduction through carbon sink, notwithstanding the large amount of their forest cover. Oita Prefecture is one of them.

The Oita Prefecture government has so far focused on energy conservation as the main mechanism for effectively reducing its CO₂ emissions [7]. Oita Prefecture’s global warming measure implementation started in 2005; it aimed at reducing energy consumption and wastes, and introducing more energy-efficient (cogeneration) systems. The result of this approach (Figure 1, [8]) does not seem to be as effective as claimed by the government. There is thus the need to reconsider this regulation and the tools used for its implementation toward the goals. Thus far, Oita Prefecture has not focused on the concept of carbon sinks for reducing its CO₂ emissions, and information is only marginally available. This is more likely because of the lack of sound scientific data needed for objective decision making about which process to use for effectively implementing the CO₂ reduction policy in the long term.

Traditional methods of estimating carbon/CO₂ sequestration through site-based methods [9,10] or using eddy covariance flux tower are known to be expensive, time consuming and limited in area coverage. Expanding the area or implementing continuous monitoring is even more costly and time consuming. Remotely sensed satellite observations have provided the scientists with an alternative method for studying the earth’s biosphere (atmosphere, vegetation, etc.). As for the vegetation, it has been demonstrated that the reflected RGB (red, green, blue), and mostly the NIR (near infrared), wavelengths contain considerable information about plants biomass [11], from which precise knowledge about the vegetation can be extracted. This gives us the possibility of biomass monitoring at lower costs and with less time loss. Today, with further research advancement, remote sensing is becoming a common analysis tool, not only for producing maps needed for categorizing the land cover type of the surface and for allocating or managing the earth’s resources [12], but also for analyzing the changes and their impacts for future land use/land cover (LU/LC) developments [13].

Using remote sensing data for evaluating regional carbon/CO₂ sequestration has been implemented in various ways, among which is combining the land cover information with the averaged carbon sequestration values of different land cover types [14]. The strength of using land cover information based on remotely sensed images is that it covers areas of regional or even global scales, making it possible to extend analysis over sites that are difficult to access on the ground. Moreover, it allows the implementation of continuous monitoring of the area, making it possible to analyze the temporal changes of the land cover.

Estimates of sequestration by forests in Japan are usually implemented at the national scale (whole Japan), focusing mostly on identifying forests as natural or planted forests. Hiroshima and Nakajima [15] focused on the potential carbon sinks in Japanese planted forests during the first commitment period of the Kyoto Protocol. Using subsidies and forest workers’ wages as predictor variables, they estimated that the privately planted forests were sequestering 8.16–8.87 MtC/yr, depending on the different scenarios of forest management and silvicultural practices. Sasaki and Kim [16] estimated both the potential and the eligible sequestrations by the forests in Japan under the condition of the Marrakesh Accords. Using the land use and the carbon stock growth models, they estimated that the forests in Japan were likely to sequester 20.1 MtC/yr (planted forests: 15.3 MtC/yr; natural forests: 4.8 MtC/yr), while under the Marrakesh Accords scheme, this was estimated to be 10.2 MtC/yr (planted forests: 7.3 MtC/year; natural forests: 2.9 MtC/yr). In general, few studies have estimated the potential of CO₂ sequestration by the forests in Japan, and only fewer considered to perform the analysis at local or regional scales, focusing more on detailed forest types. This can be easily understood: the Kyoto Protocol deals with national level estimates, so many studies are involved in the estimates at that scale.

However, we have to note that the estimates by those authors mentioned above are based on modeling, while the information of where and how much of those resources exist is not considered. This is a serious weakness, especially when forest resources management is concerned. Moreover, we also have to consider that forests are more diverse than just being natural or planted forests. Because both tree types and tree ages make differences in CO₂ sequestration [9,17,18], and for a more realistic estimation of the sequestration, we emphasize that forests have to be analyzed on a local or regional scale in order to get more precise and detailed classes for each forest category (tree type) as well as for each tree age range. Based on this detailed information about the forests, precise estimates of the forests sink capacity can be obtained. This is also a way of providing reliable information for sound decision making for sustainable forest management and development.

The objective of this research is to estimate CO₂ sequestration by forests type in Oita Prefecture, based on precise evaluation of forests extent in Oita. This will make it possible to quantitatively estimate their potential storage capacity, and therefore maximize the usage of this environmental resource as a CO₂ sink. These findings will help the Oita Prefecture g to devise better scenarios for mitigating global warming. Local and regional scale analysis of the land cover and estimation of CO₂ sequestration will eventually contribute to the updated information of Japan’s CO₂ reduction target achievement.

The study site, datasets and methodology are presented in Sections 2–4, respectively. Geographic Information Systems (GIS) and remote sensing technology are used for the analysis. The analysis results are given in Section 5, followed by discussion in section 6 and concluding remarks in Section 7.

2. Study Site

The research area, Oita Prefecture, is located on the Kyushu Island in the southwestern part of Japan (Figure 2), between 130.81°E and 132.11°E, 32.72°N and 33.75°N, with a total land area of 6,339 km². Official figures show that 4,540 km² (or 72%) of Oita Prefecture’s lands are covered by forests [19].

Oita Prefecture is also known for its rich and diversified Onsen (hot springs), found especially in Beppu City. A volcanic belt (the Kirishima Volcanic Range), with several mountains ranges, culminating up to 1,758 m in altitude, runs from north to south through the prefecture, contributing to the existence of the many hot springs.

Climatic conditions divide the prefecture in two main areas: (1) the area affected by the monsoon winds coming from the coast of the Japan Sea, and bringing in large amounts of rainfall and large number of rainy days in winter, and (2) the area characterized by large amounts of rains in summer resulting from the moisture flowing in from the Pacific Ocean, and by dry and fine weathers in winter. Because of these wet climate conditions, annual rainfall in Oita Prefecture ranges from 1,500–1,600 mm over coastal areas, to 2,500–3,000 mm in the mountainous areas.

3. Data Sets

Two optical satellite data sets from the Land Remote Sensing Satellite Enhanced Thematic Mapper, Plus (LANDSAT ETM+) with 30 m spatial resolution published by the United States Geological Survey (USGS) and two microwave satellite data sets from the Advanced Land Observing Satellite (ALOS) Phased Array L-band Synthetic Aperture Radar (PALSAR) provided by Japan Aerospace Exploration Agency (JAXA), were selected for our analytical objectives. The selection of the data sets was done by considering the following three critical issues. First, we had to choose a cloud-free image for analysis. In fact, clouds block the solar radiation from reaching the ground, rendering the covered area unavailable for the optical observing platform. Second, the Scan Line Corrector (SLC) failure that occurred after 31 May 2003 makes the LANDSAT ETM+ line of sight trace a zigzag pattern along the satellite ground track. As a result, a gap fill has to be performed in order to correct the data so as to become usable for analysis [20]. Third, we had to select an image old enough for us to create a digital LU/LC map (or digital vegetation model, here after DVM) of Oita Prefecture that could be used later on as a reference digital map for further quantitative LU/LC change analysis in this area. As a result, we selected the cloud-free LANDSAT ETM+ images captured on 25 May 2002, and 30 December 2006.

For the microwave data, we chose the level 1.5 product of ALOS PALSAR L-band data of 15 September 2009 and 14 October 2009, to cover the study area with an ascending Fine Beam Dual (FBD) Polarization, characterized by 34.3 degree of off nadir angle and 12.5 m of ground range pixel spacing, as the best fit for our analysis purposes.

The National Survey on the Natural Environment [21], published by the Ministry of Environment of Japan, provides the status of lands, surface waters and coastal areas all over the nation. The survey results are compiled and published in the form of written reports, downloadable maps, etc. The vegetation map so provided is utilized for identifying specific LU/LC sites for ground truth training data needed for land features classification and for the accuracy assessment of the classified images. Field observations were also performed for collecting additional ground truth data and for the identifying land features for the selected sites all over the prefecture.

For the correction of the remote sensing data, we utilized the digital model of terrain surface published by the Earth Remote Sensing Data Analysis Center (ERSDAC) available from the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) Global Digital Elevation Model (GDEM) with a spatial resolution of 30 m and an elevation accuracy of ±7 m. The data were downloaded from the websites and concatenated so as to cover the study area and to make them ready for further use in the analysis.

4. Analysis Method

First, gap fill was performed before the correction of the images for the LANDSAT data using the NASA Goddard Space Flight Center frame and fill tool [22]. Then, the Integrated Radiometric Correction (IRC) method developed by Kobayashi and Sanga-Ngoie [23] was applied to the selected LANDSAT ETM+ images in order to correct the radiometric distortions caused by both the topography and the atmosphere, so as to obtain the actual (corrected) irradiance value. The IRC method is a comprehensive radiometric correction method aimed at simultaneous correction of atmospheric, solar, and topographic effects inherent to remote sensing data. The IRC makes it possible to retrieve the at-surface reflectance and radiance values as if the atmosphere were totally transparent and the underlying surface absolutely flat, using

$$L_{\text{g}}=(L_{\text{s}}^{*}-L_{\text{h}})\frac{A}{0.5(1+T_{\text{r}})T_{\text{r}}{T}_{\text{w}}^2}$$

(1)

where, L_g is the corrected radiance from the surface; L_s^* is the at-satellite observed spectral radiance from a sloped terrain and L_h is the upwelling atmospheric spectral radiance also known as haze; A is the topographic correction factor; T_r and T_w are atmospheric transmittance functions: Rayleigh scattering transmittance function (T_r) and water vapor transmittance function (T_w). The Dark Object Subtraction Method (DOSM) was used for the removal of haze (L_h). Haze values were computed for each band using the reflectance values over deep water areas selected as areas of known zero reflectance in the image [24].

Table 1. Computed haze (L_h) values for each band using the Dark Object Subtraction Method (measured in Wm⁻²·μm⁻¹·sr⁻¹).

	Band 1	Band 2	Band 3	Band 4	Band 5	Band 7
May	40	11.8	5	1.3	0.13	0.05
December	26	11.8	5	1.3	0.13	0.04

shows the computed haze values for each band. The decrease of haze values with increasing wavelengths is remarkable.

A hybrid classification method [25,26] is used here in the classification process for producing the digital LU/LC map of Oita. Sanga-Ngoie [27] stressed on the need to produce a detailed digital LU/LC map or DVM, especially for the vegetation categories, in order to obtain more precise estimations of CO₂ sequestration values. Here, the May data will be the main data set for the analysis while the December data and the PALSAR data will be supplemental data for additional information to the classification/reclassification process so as to produce more detailed land cover maps.

For more precise estimation of CO₂ sequestration, the PALSAR data were used for estimating stem volume, from which the information about tree age could be deduced. The PALSAR Digital Number (DN) values were converted to the backscattering intensity (σ⁰) in decibel units using Equation (2),

$${\sigma }^0=10{\text{log}}_{10}{\text{DN}}^2+\text{CF}$$

(2)

where CF is the calibration factor for the PALSAR products, for the level 1.5 product CF = −83.0 dB [28]. For a correct interpretation of the backscatter signatures, corrections for the effect of local incidence angle and normalization for the true pixel area were necessary. The corrected backscatter in gamma nought (γ⁰) format can be obtained from the sigma nought (σ⁰) value (Equation (3)) according to Castel et al.[29] and Ulander [30]. This approach has also been applied by other authors [31–34].

$${\gamma }^0={\sigma }^0\frac{A_{\mathit{flat}}}{A_{\mathit{slope}}}{\left(\frac{\text{cos}{\theta }_{\mathit{ref}}}{\text{cos}{\theta }_{\mathit{loc}}}\right)}^n$$

(3)

where θ_loc and θ_ref represent the local incidence angle and a reference angle for the normalization of the backscatter (e.g., off nadir angle), respectively. A_slope and A_flat represent the true pixel area and the local pixel area for a theoretically flat terrain, respectively. For bare surfaces, the exponent n is equal to 1. For vegetated surfaces, n expresses the variation of the scattering mechanism due to the presence of a volume on the sloped terrain. It is thus being related to the optical depth of the vegetation, and is estimated to be lower than or equal to 0.36 for the Japanese Earth Resources Satellite-1 (JERS-1) L-band SAR backscatter. We opted to use the 0.36 value, because the scattering properties of the forests should be taken into account. For future works, more precise estimates should be reconsidered, especially in the case of different polarizations.

The overall flowchart of the methodology is given in Figure 3, and consists of four main steps described hereafter. The Clark Lab’s IDRISI Taiga GIS platform [24] was used for all the analysis.

First, unsupervised classification [24,25,35] was performed to obtain a cluster image to be used for identifying the ground cover types. The IDRISI Taiga ISOCLUST module was used for the clustering process, yielding 15 distinct clusters.

Training sites were then constructed upon the cluster image and used for the supervised classification process based on the Maximum Likelihood algorithm. Data from the National Survey on the Natural Environment [21], together with our own field investigation data, were used as ground truth references for land features identification. The strength of this method is that the need of delineating unimodal spectral classes beforehand can be overcome based upon a representative subset of image data [25]. The IDRISI Taiga MAXLIKE (maximum likelihood) module and the LANDSAT ETM+ 12 bands (bands 1–5 and 7 (May and December)) along with Normalized Difference Vegetation Index (NDVI) data of May and December were used in order to perform the supervised classification. This yielded a classified image with 12 categories in total. These were reclassed later on to give our final product: a classified image with 10 LU/LC categories, designed with a special attention to the forests types found over the analysis area.

In the second step, accuracy assessment was made, based on the final image. Stratified random sampling of the image yielded 29,000 sample points used for the accuracy assessment. Using the ERRMAT module of IDRISI Taiga, the error matrix was built in reference to the land features on the vegetation map from The National Survey on the Natural Environment data [21] as ground truth.

In the third step, the information obtained from the reclassed image was used together with the CO₂ sequestration rates per unit area [9] for the estimation of the sequestration values by each forest type. For the coniferous and deciduous broadleaf forests, we performed the estimations not only by forest type but also by breaking these forest categories down into categories per tree age in the final steps of the work.

The classification of the coniferous forests per tree age was done by calculating the NDVI values using Equation (4), where IR stands for the infrared band (band 4) and R for the red band (band 3). The relationship between NDVI and tree age was deduced by Ishii [36] based on the Forest Tree Age growth curve and multitemporal NDVI data. Referring to this growth curve, we could thus classify the coniferous forests area into separate tree age categories.

$$\text{NDVI}=(IR-R)/(IR+R)$$

(4)

Furthermore, considering the fact that NDVI saturates with age, especially for trees older than 15 years, we attempted to improve the classification of forests area (especially coniferous and deciduous broadleaf forests) into separate tree age categories, using the backscattering intensity information derived from the PALSAR data. This was implemented for comparing the results with the tree age estimation based on the NDVI method. The relationship between backscattering intensity and stem volume, or biomass, has been analyzed by various researchers [31,34,37–42]. These various studies state that the backscattering signature correlates with the forest parameters, and that SAR images can thus be used to extract information related to the forest biomass. Because we did not collect ground data of the stem volume in our study area, we opted to use the polynomial model developed by Wijaya [40] (Equation (5)) to estimate the stem volume (m³/ha).

$$V=256.85-(65.458{\sigma }^0)-(6.8{\sigma }^{0^2})$$

(5)

where, V is the stem volume (m³/ha) and σ is the backscattering intensity. Notwithstanding the fact that this expression was developed for forests in Indonesia, we decided to apply this model because Wijaya [40] used it for very similar data to ours (PALSAR L-band FBD 34.3°HH/HV specifications). The stem volume map produced using this model was then classified into ages among the forests, referring to the relationship of the tree age and stem volume information provided by the Oita Prefectural Government. The classification was done for the coniferous and deciduous broadleaf forests, but not for the evergreen broadleaf forests for which no stem volume data were available. Sugi (cedar) tree data and Kunugi (saw tooth oak) tree data were used as representative for the coniferous forests and deciduous broadleaf forests, respectively.

5. Results

5.1. Radiometric Correction

Quantitative evaluation of the performance of the correction is done by assessing the absence of correlation between cos i (i: solar incidence angle), which is a measure of the solar illumination on the sloped area, and the corrected radiance data. The regression analysis between cos i and noncorrected data (Figures 4 and 5(1)) shows that these parameters are highly correlated: higher radiances for higher cos i values, because of the various inherent errors, including the atmospheric and topographic effects. However, lower correlation, or mostly no correlation, between those two parameters is obtained for radiometrically corrected radiances (Figures 4 and 5 (2)). Figure 6 shows an example of the difference of images before and after correction, respectively, for the band 4 image of May data. The zoomed-in area is a mountainous location, showing that the radiance values are corrected in these mountainous areas, while almost no changes are detected over the flat residential areas on the top left side of the inset.

The use of radiometrically corrected data ensures that more precise land use/land cover classification of optical satellites data will be performed, with minimum perturbations related to atmospheric effects (haze, transmittance, etc.) or topographic effects (hillshading).

5.2. Oita LU/LC Digital Map

The LANDSAT ETM+ image captured on the date of 25 May 2002 was classified referring to the Natural Survey on the Natural Environment vegetation map [21] and our own field investigation records as ground truth, into the following 12 land cover types (

Table 2. Land cover classes: (a) Original outputs; (b) The final product.

): cedar–cypress (Cryptomeria japonica–Chamaecyparis obtusa), saw tooth oak (Quercus acutissima), live oak (Castanopsis Cyclobalanopsis), silver grass, grasslands, bare lands, water bodies, crop lands, paddy fields, residential areas, deciduous broadleaf forests and evergreen broadleaf forests. These detailed classes were then reclassed into the following final 10 LU/LC categories (Table 2(b)): coniferous forests, silver grass, grass lands, bare lands, water bodies, crop lands, paddy fields, residential areas, deciduous broadleaf forests and evergreen broadleaf forests (Figure 7).

Table 3. The error matrix for the classified LU/LC map. Area for ground truth categories, together with the Commission Errors (ErrCom), are given in the columns, while the areas for the classified categories and the Omission Errors (ErrOm) are given in the rows.

Reference
		Coniferous	Silver Grass	Grasslands	Bare lands	Water bodies	Crop lands	Paddy fields	Residential	Deciduous-Bl	Evergreen-Bl	Total	ErrCom	Users A
Classified	Coniferous	9,305	100	52	16	61	125	158	39	1,086	1,543	12,485	0.2547	74.53%
	Silver Grass	6	365	48	7	2	7	3	1	17	0	456	0.1996	80.04%
	Grasslands	69	192	211	17	8	86	71	26	65	18	763	0.7235	27.65%
	Bare lands	5	10	14	41	12	5	7	59	3	2	158	0.7405	25.95%
	Water bodies	5	9	2	12	279	0	8	6	1	14	336	0.1696	83.04%
	Crop lands	179	29	67	52	28	235	237	197	42	71	1,137	0.7933	20.67%
	Paddy fields	309	128	129	109	223	216	1,339	349	120	105	3,027	0.5576	44.24%
	Residential	5	15	40	55	70	34	104	651	7	9	990	0.3424	65.76%
	Deciduous-Bl	1,190	363	129	18	4	85	73	12	1,788	500	4,162	0.5704	42.96%
	Evergreen-Bl	1,734	77	89	25	16	257	108	44	1,123	1,729	5,202	0.6676	33.24%
	Total	12,807	1,288	781	352	703	1,050	2,108	1,384	4,252	3,991	28,716
	ErrOm	0.2734	0.7166	0.7298	0.8835	0.6031	0.7762	0.3648	0.5296	0.5795	0.5668		0.4448	Overall Accuracy
	Producers A	72.66%	28.34%	27.02%	11.65%	39.69%	22.38%	63.52%	47.04%	42.05%	43.32%			15,943/28,716 = 55.52%

shows the results of the error matrix of the accuracy assessment obtained by using the ERRMAT module of the IDRISI Taiga GIS software, while

Table 4. Fuzzy Error Matrix [43,44] following the same methodology as the traditional error Non-diagonal cells in the matrix contains accepted and poor columns. The numbers in the accepted estimating fuzzy accuracy, the numbers in the poor columns are considered as an error. The fuzzy percentage of sites where the acceptable reference labels matched the classified label.

Reference																					Users Accuracy
	Coniferous		Silver Grass		Grasslands		Bare lands		Water bodies		Croplands		Paddy fields		Residential		DBL		EBL		Total Fuzzy	Percent Fuzzy
	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor	Accept	Poor
Coniferous	9,305	0	0	100	0	52	0	16	0	61	0	125	0	158	0	39	642	444	911	632	10,858/12,485	86.97%
Silver Grass	0	6	365	0	27	21	0	7	0	2	4	3	0	3	0	1	0	17	0	0	396/456	86.78%
Grasslands	0	69	108	84	211	0	0	17	0	8	47	39	39	32	0	26	0	65	0	18	405/763	53.1%
Bare lands	0	5	0	10	0	14	41	0	0	12	0	5	4	3	25	34	0	3	0	2	70/158	44.07%
Water Bodies	0	5	0	9	0	2	0	12	279	0	0	0	3	5	0	6	0	1	0	14	282/336	84.07%
Croplands	0	179	22	7	43	24	0	52	0	28	235	0	199	38	182	15	0	42	0	71	681/1,137	59.92%
Paddy fileds	0	309	99	29	98	31	79	30	140	83	185	31	1,339	0	290	59	0	120	0	105	2,229/3,027	73.62%
Residential	0	5	0	15	0	40	23	32	0	70	22	12	66	38	651	0	0	7	0	9	762/990	76.97%
DBL	934	256	0	363	0	129	0	18	0	4	0	85	0	73	0	12	1,788	0	443	57	3,165/4,162	76.04%
EBL	1,194	540	0	77	0	89	0	25	0	16	0	257	0	108	0	44	999	124	1,729	0	3,922/5,202	75.39%
Producers	11,433/12,807		594/1,288		379/781		143/352		419/703		493/1,050		1,650/2,108		1,148/1,384		3,429/4,252		3,083/3,991		Fuzzy Overall Accuracy
Accuracy	89.27%		46.12%		48.53%		40.63%		59.60%		46.95%		78.27%		82.95%		80.64%		77.25%		22,771/28,716 = 79.30%

shows the fuzzy error matrix modified from the conventional error matrix based on the concept of fuzzy accuracy assessment developed by Green and Congalton [43,44]. The error matrix so obtained tabulates the errors of omission and commission, as well as the marginal and total errors. The overall accuracy of the produced LU/LC map is 55.5% and the fuzzy overall accuracy is 79.3%.

5.3. CO₂ Sequestration by Forest Types

The three dominant forest types found in Oita Prefecture are: coniferous, deciduous broadleafs and evergreen broadleaf forests. Estimations of CO₂ sequestrations per forest type are obtained by multiplying the observational CO₂ sequestration rates [9,27,45] by the total area of each forest type (

Table 5. CO₂ sequestration by forest types in Oita Prefecture using the IPCC Method.

Forest Type	CO2 Sequestration (tCO2/ha/yr) [9,47]	Area (ha)	Total CO2 Sequestration (tCO2/yr)	Total C Sequestration (tC/yr)
Coniferous	13.5 ± 4.2	264,152	3,566,052	972,560
Deciduous-BL	8.7 ± 3.0	89,208	776,109	211,666
Evergreen-BL	18.1 ± 3.0	124,494	2,253,341	614,548
Total		477,854	6,595,502	1,798,774

). The following results are obtained for the total CO₂ sequestrations per forest type: 3.57 MtCO₂/year for the coniferous forests, 0.78 MtCO₂/year for the deciduous broadleaf forests, and 2.25 MtCO₂/year for the evergreen broadleaf forests. This sums up into a total CO₂ sequestration of 6.60 MtCO₂/year for all the forested areas in Oita Prefecture (Table 5). These findings emphasize the fact that more precise sequestration estimates can be obtained by using detailed LU/LC maps (with detailed forests categories and precise areas) than by just dividing the forest types into natural and planted forests only as recommended by the IPCC [46].

Carbon sequestration equivalents are simply obtained by multiplying the CO₂ sequestrations by 12/44, where 12 and 44 are the molecular masses of the carbon and the CO₂, respectively. Our aim is to evaluate atmospheric CO₂ sequestration by forests. We also present in the tables the carbon sequestration equivalents, for informational purposes.

5.4. CO₂ Sequestration by Tree Ages (NDVI Method)

Further estimation of the CO₂ sequestration was considered for the coniferous forests. These forests were classified into four different categories (Figure 8), based on the NDVI values [36], each category indicating the range of the tree age for the coniferous forests: 0–5, 5–10, 10–15 and over 15 years old categories. Based on this detailed classification, a total CO₂ sequestration of 2.35 MtCO₂/yr was obtained for coniferous forests in Oita (

Table 6. CO₂ sequestration by coniferous forests (cedar–cypress) per tree age using the NDVI Method.

Coniferous Age	CO2 Sequestration (tCO2/ha/yr) [47]	Area (ha)	Total CO2 Sequestration(tCO2/yr)	Total C Sequestration (tC/yr)
0–5	0.6	15	9	2
5–10	4.8	1,316	6,317	1,723
10–15	13.5	11,772	158,922	43,342
15–50	8.7	251,048	2,184,118	595,669
Total		264,151	2,349,366	640,736

).

The averaged sequestration value used for planted forests following the traditional (IPCC) methodology [46] is 13.5 tCO₂/ha/yr. We have to note that in Table 6 (sequestration value per unit area per tree age), a similar value of 13.5 tCO₂/ha/yr is shown to be valid only for trees 10–15 years old. This age range is also the period of time when the coniferous trees seem to have the highest growth rate. There is therefore a strong possibility of overestimating the sequestration results if the estimates are made using these averaged values (which turn out to be the values for the tree age category with the highest sequestration rate) for all the tree age categories. In fact, the comparison of the total sequestration values obtained for coniferous forests by using the averaged value (3.57 MtCO₂/yr) with our results in Table 6 (2.35 MtCO₂/yr) shows a 34.2% overestimation in CO₂ sequestration for the former methodology.

5.5. CO₂ Sequestration by Tree Age using PALSAR Data (Stem Volume Method)

Using Equation (5) and the backscattering intensity derived from the PALSAR data, we calculated stem volumes for forests in Oita, as shown in the stem volume map (Figure 9). The stem volume values for various tree ages as provided by the Oita Prefectural Government are shown by the graphs in Figure 10. Then, the so-obtained stem volume map (Figure 9) is reclassified into range of ages, by considering the different growth rates among the various regions of the study area, based on Figure 10.

Coniferous forests data were analyzed and the results compared with those obtained by using the NDVI Method. The same was done for the deciduous broadleaf forests. Here, obtained forest ages are more detailed or more specific than those used for the NDVI method, because of the difference between the saturation points for the NDVI and the stem volume methods. In fact, the NDVI values saturate after 0.7, showing similar values for all the ages above, therefore the difficulty in classifying trees of more than 15 years of age. For the stem volume method, the same can be said except that here, the saturation point is higher than for the NDVI, as is widely discussed in the literature [40,47].

Table 7. CO₂ sequestration of coniferous and deciduous broadleaf forests per tree age using the PALSAR data and the stem volume method.

Coniferous Age	CO2 Sequestration (tCO2/ha/yr) [47]	Area (ha)	Total CO2 Sequestration (tCO2/yr)	Total C Sequestration (tC/yr)
0–5	2.9	3,184	7,323	1,997
5–10
10–15	13.5	8,941	112,657	30,725
15–20	14.8	90,441	1,293,306	352,720
20–25	11.3	147,084	1,500,257	409,161
25–30	10.2	5,358	46,615	12,713
30–50	6.4	662	3,972	1,083
Total	11.59	255,640	2,964,130	808,399
Deciduous-Bl Age
0–5	5.5	310	1,705	465
5–10	9.4	614	5,772	1,574
10–15	8.2	1,205	9,881	2,695
15–20	8.2	1,557	12,767	3,482
20–25	5.9	2,168	12,791	3,488
25–30	5.8	2,709	15,712	4,285
30–50	3.2	78,044	249,741	68,111
Total	3.56	86,607	308,369	84,100

shows the new total CO₂ sequestration estimates for coniferous and broadleaf forests per age category, respectively. In general, the new total CO₂ sequestration for coniferous (2.96 MtCO₂/yr) is found to be larger than that obtained by using the NDVI method (2.35 MtCO₂/yr) (Table 6). As for the deciduous broadleaf forests, a total sequestration of 0.38 MtCO₂/ha/year is obtained (Table 7) in comparison with the 0.78 MtCO₂/ha/yr (Table 5) by using the traditional averaged method. These differences come from the fact that the sequestration values, which were averaged for all the trees of the same category (traditional averaged method), or for all the trees of more than 15 years of age due to NDVI saturation (NDVI method), are now broken down into detailed age categories up to the year 50 (stem volume method).

In fact, by using the detailed sequestration rates per tree ages (stem volume method), sequestration overestimates that implicitly occur with the averaged method [46] for coniferous trees of less than 10 years of age and for trees of more than 20 years of age, as well as the underestimation for 15–20 years old trees, can altogether be avoided. The same can also be said, mutatis mutandis, for deciduous broadleaf forests for which the 5–10 years category is shown to have the highest sequestration rate (9.4 tCO₂/ha/yr) according to the stem volume method (Table 7), as compared to the averaged values for all ages (8.7 tCO₂/ha/yr) in Table 5. In general, the highest sequestration rates are shown for 10–25-year-old coniferous forests, and for 5–20-year-old deciduous broadleaf forests.

6. Discussion

6.1. Land Cover Map (DVM)

Focus has been put on producing a land cover map in which land cover types are defined not just simply as forested and nonforested areas [46], but more in detailed forest subclasses and other nonforested areas. Sanga-Ngoie [27] had stressed the need to produce a detailed digital land cover map, especially for the vegetation categories, in order to allow for more precise estimations of CO₂ sequestration values.

Various combinations of approaches and data sets (noncorrected vs. corrected, single vs. multitemporal) were examined in order to produce the best results for a highly accurate land cover map (Figure 7). Difficulty did arise when using single temporal datasets for producing land cover maps, especially when one has to classify the vegetation types in detail. For example, the spectral signatures (mean radiance values per training site for each band) for forested land cover types given in Figure 11 show that it is very hard to differentiate the evergreen broadleaf forests and deciduous broadleaf forests categories in the May data, because of similar signatures for most of the bands in this period of the year. However, it is shown that clear discrimination among all the different forest types can be obtained from the December data (bands W1–W5 and W7). In general, although the radiance values in December (winter season) are lower than those in May (spring season) for all the forests types, and coniferous clearly differentiated from other broadleaf forests in both seasons, it can be said that, when seasonal changes in the forests canopy are considered: (1) clear differentiation among the deciduous and evergreen subgroups of broadleaf forests can be obtained with larger radiance values for deciduous broadleaf forests than for evergreen broadleaf forests, (2) discrimination within each broadleaf forest subgroup can be obtained, with radiance values (in Wm⁻²·μm⁻¹·sr⁻¹) for saw tooth oak and live oak forests larger than those for other deciduous and evergreen broadleaf forests, respectively. Similar improvements are also established in terms of accuracy assessment, showing by far larger accuracy values (72.3%, 42.1% and 43.3%) in Table 3 for multitemporal corrected images than the 40%, 14% and 22% obtained for coniferous, deciduous and evergreen broadleaf forests, respectively, based on the analysis of single and noncorrected bands. These findings firmly establish the importance of considering seasonal differences (multitemporal data sets) in order to improve the land cover map results, because they make it possible to distinguish various land cover categories normally difficult to classify from one single temporal data sets.

6.2. Accuracy Assessment

It is a common practice in the remote sensing society to use a conventional error matrix for checking the reliability of the classified land cover map. Detailed information from within the error matrix is often overlooked because people mostly focus only on the overall accuracy. Congalton and Green [48] recommended to consider the error matrix as a jumping-off point for identifying not simply the error itself, but mostly the sources of confusion that might come from the differences between the map created from remotely sensed data and the reference data, or from erroneous human interpretation, for example when defining these features from aerial photos. The confusion might thus lead to errors in the process of producing reference maps, which would in their turn produce errors when used later on as a ground truth references. Accuracy assessment should thus be a process for identifying errors and understanding the problem of the misclassification so as to prevent them. We stress that more attention has to be paid to this issue instead of concentrating on the overall accuracy alone.

In this study, fuzzy accuracy assessment [43,44,48] was applied along with the conventional accuracy assessment process in order to compare how the accuracy varies when we fuzz the boundaries of the classes. Vegetation Condition Index (VCI) [49] was used to compute the degree of acceptance for each class. Considering the complexity of the real environment, we suggest that fuzzy accuracy assessment should also be implemented in the same way as the conventional accuracy assessment for any sound classification scheme.

6.3. NDVI, PALSAR data and CO₂ Sequestration Estimates

The overall sequestration results have been summarized in

Table 8. Summary of CO₂ sequestration values and total sequestrations for the different methods of estimation.

Estimation Method	Coniferous (tCO2/ha/year)		Deciduous-BL (tCO2/ha/year)		Evergreen-BL (tCO2/ha/year)	Total CO2 Sequestration (tCO2/year)	Total C Sequestration (tC/year)
Conventional	13.5		8.7		18.1	6,595,502	1,798,774
NDVI	Age 0–5	0.6	8.7		18.1	5,378,816	1,466,950
	5–10	4.8
	10–15	13.5
	15–50	8.7
Stem Volume	Age 0–5	2.9	Age 0–5	5.5	18.1	5,525,840	1,507,047
	5–10		5–10	9.4
	10–15	13.5	10–15	8.2
	15–20	14.8	15–20	8.2
	20–25	11.3	20–25	5.9
	25–30	10.2	25–30	5.8
	30–50	6.4	30–50	3.2

so as to enhance our findings in this paper. As can be seen in this table, if the conventional IPCC method is used alone, coniferous and deciduous broadleaf forest sequestrations are set to be 13.5 tCO₂/ha/yr and 8.7 tCO₂/ha/yr, respectively. The NDVI methodology made it possible to break the Coniferous sequestration down to each age category, showing very clearly the range of overestimation implied when the IPCC method (with no special attention to tree age) is used, especially for trees younger than, or older than, 10–15 years of age. However, because of the saturation effect for the NDVI values, sequestration rates could not be discriminated for trees older than 15 years with the NDVI method.

The use of PALSAR data made this discrimination possible, yielding sequestration rates per tree age, of up to trees of 50 years of age, not only for the coniferous, but also for the Deciduous broadleaf forests. It is thus clearly established that sequestration rates of 13.5 tCO₂/ha/yr (cfr IPCC Method) are valid only for 10–15-year-old coniferous forests, while deciduous broadleaf forests reach, or go beyond, the 8.7 tCO₂/ha/yr threshold only during the 5–10-year-old range. Any use of these IPCC rates for trees of different ages from the classes mentioned here above will lead to sequestration overestimation. The same is also true for the results obtained using the NDVI method for the coniferous, except for the underestimation that occurs for the trees between the 15 and 20-year range.

It is very important to note that, the averaged sequestration value for the coniferous forest type matches the sequestration value per tree ages obtained here for the 10–15-year-old coniferous trees. This is also the category of the most fast-growing period for coniferous trees. The same is also true for the deciduous broadleaf forests for the 5–10-year-old trees. This strongly implies that, if estimation is made based on the averaged values alone, it will round up the values for other trees in different growth stages, giving way to an overall overestimation of the carbon sink capacity, both for natural or planted forests [46]. Our results emphasize the fact that sequestration differs among tree ages: this has to be considered, otherwise erroneous values of forests’ sequestration capacity will be obtained, and a number of negative consequences in the further use of these values can be foreseen.

This demonstrates and emphasizes the importance of using PALSAR data (stem volume method), based on the estimation of the stem volume computed from the PALSAR backscattering intensity data by applying the model developed by Wijaya [40] for Indonesian forests. Nowadays, many works attempt to analyze the correlation between backscatter with different polarizations and the stem volume or biomass for forests with the same tree types [41,42]. They do not look through the backscattering signatures among different vegetation types for different polarizations. It is known that the backscattering signature gets affected by the soil moisture or the topography and it can also be argued on how well Wijaya’s model applies to the study area in Oita because of different climate conditions and vegetation. The impacts, or the eventual limitations, on the applicability of Wijaya’s expression to forests in Oita will be one of our most urgent tasks to deal with in the future.

6.4. CO₂ Sequestration Estimates by Oita Government

Recently, the Oita Prefectural Government published its sequestration estimates by forests of Oita [8]. They divided Oita forests into natural and planted forests, and calculated the sequestration values as prescribed in the IPCC Guideline [46], using the following equation:

$$\begin{array}{l}\text{Annual}\hspace{0.17em}\text{Carbon}\hspace{0.17em}\text{Sequestration}\hspace{0.17em}(\text{per}\hspace{0.17em}\text{ha})=\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\text{Annual}\hspace{0.17em}\text{Growth}\hspace{0.17em}({\text{m}}^3/\text{ha})\times \text{Magnification}\hspace{0.17em}\text{Factor}\times \text{Bulk}\hspace{0.17em}\text{Density}\times \text{Carbon}\hspace{0.17em}\text{Content}\end{array}$$

(6)

where the magnification factor is a factor for expanding the stem volume into an estimate of the whole volume of the tree, including the volume of all leaves, branches and roots together. A magnification factor of 1.7 and 1.8 is proposed for coniferous trees and broadleaf trees, respectively. Bulk Density represents the weight per cubic meter of volume. We have to note that, although a common value is used everywhere in Japan, we presume that the values of these two factors should vary among trees according to, not only their types, but also their ages and the regions of growth. The carbon content is usually set to be 0.5 [46], meaning that the composition of wood contains roughly 50% of carbon. CO₂ sequestration is simply obtained by multiplying the C sequestrations by 44/12, where 12 and 44 are the molecular masses of the carbon and the CO₂, respectively. Finally, in order to obtain the total carbon sequestration by the forests, Equation (6) will have to be multiplied by the forested total area (in ha).

Published Oita Prefecture’s CO₂ sequestration results [8] using the IPCC Method [46] are still not very conclusive and need some more work and systematic streamlining. Strong discrepancies are found among the few results published so far, more likely due to both the difference in forest classification and in the sequestration assumptions used for the calculations. We have to stress here that these prefectural values are obtained without consideration to the differences in sequestration rates by tree type or tree age, they are also relatively low compared to the values in the literature [9,27,45], and to our present results (Tables 6–8) in which we emphasized the diversity in sequestrations among tree type and age.

In summary, it is shown that Oita Prefectural Government is underestimating its forests sequestration capacity. We would like to recommend the prefecture to consider the methodology developed here in the case that the carbon sink concept is to be considered as part of the prefectural CO₂ countermeasure plan; because the government’s present estimates could lead to incorrect decision making which ought to be avoided.

6.5. Decision Making in Forests Management

Our results have shed new light on the role and the potentials of forests in Oita Prefecture concerning CO₂ sequestration. If well implemented, the present methodology is more likely to help the prefecture in devising better scenarios for forest resource management, and for global warming mitigation, as well. Through this methodology, resources could be localized and their areal covers quantitatively estimated using the satellite optical remote sensing data, while the status of the resources could be analyzed using the SAR data.

By comparing with the observed sequestration rates per tree age, it was shown that the highest sequestration rates are obtained for 10–25-year-old trees for coniferous, and for 5–20-year-old for deciduous broadleaf forests. Tree harvesting time and afforestation plans can thus be pinpointed for specific forested areas, taking into consideration the important information about sequestration rates per tree type and per tree age. Planted trees will thus need to be checked if they can yield good timber products at around 20–25 years or 15–20 years of age for coniferous or for deciduous broadleaf forests, respectively. Bibilis et al.[50] emphasize the importance of having more precise information about timber products, because the percentage of acceptance to grade requirements is known to increase with tree age. Altogether, these findings can be a strong information source in support of sound decision making and actions not only for the Oita Prefectural Government, but also for other prefectures in Japan, or all over the world.

7. Concluding Remarks

This research has employed multitemporal LANDSAT ETM+ and ALOS PALSAR data together with field data collected in Oita Prefecture, Japan, to estimate the CO₂ sequestration potential of the forests in Oita. We have shown that sequestrations obtained by using averaged values for natural or planted forests without forest-type and tree-age discrimination can lead to large overestimations of the sink capacity of forests. It has been made clear here that, by combining the optical and SAR remote sensing data, it became possible to obtain more realistic CO₂ sequestration amounts by the forests, based on detailed discrimination of forest type and tree age. Moreover, the detailed information about the spatial location and the status of the forest resources is more likely to provide the decision makers with a sound basis for more effective policies and actions aiming at the sustainable use and conservation of the environmental resources.

In this study, the retrieval of tree ages using the stem volume information extracted from SAR data was implemented only for coniferous and deciduous broadleaf forests, taking cedar (Sugi) and saw tooth oak (Kunugi) trees as representatives of those forest types, respectively. With this method, we have obtained more realistic sequestration amounts than with the other methods. This methodology can thus help to investigate more in-depth the CO₂ sequestration values for trees of different ages and types, especially for the coniferous or other evergreen or deciduous broadleaf forests, which are found over the study area. We foresee that this would lead to even more realistic estimation of the CO₂ sequestration of the real environment.

The correct positioning of PALSAR images is another general and important issue, especially when these images are to be fused with optical data. This problem can lead to errors, especially when constructing the training sites. Because of the positioning problem, the training sites coordinates in the optical and the SAR data will not match, corrupting thus the spectral signatures information for each land class. Orthorectification is usually used to deal with this issue, but it was not implemented in this research. The correction of the PALSAR image using the methods by Castel et al.[29] does not account for layover and shadowing effects. This leads to some basic problems with the SAR data, especially in Japan, a country with many mountains and slopes. Small [51] has developed a new method for improving the normalization of the SAR image, in which the layover and shadowing effects are considered. We intend to include this method in our future work in order to overcome these issues and improve the results obtained in the present work. In so doing, we will have room for more options for the satellite data to be acquired, such as using images for ascending or descending modes, or for any other modes available.