Modeling of annual net primary production of a forest in the Taramakau Valley... - natural-foundation-science.org

1. Introduction

Over the past 30 years, a considerable number of global vegetation models have been developed to investigate the land carbon sink, which is one of the major concerns on the implementation of Kyoto Protocol. Usually, an essential part of these models is the sub-model of simulating net primary production (NPP). Net primary production is the mainly (even only) source of carbon sink in the terrestrial ecosystem, which explains most annual carbon fluxes between the atmosphere and biosphere. Most NPP simulation models are based on the interactions between vegetation growth and its available resources, such as radiation, available soil water, soil nutrient, and ambient temperature.

Among these resources, radiation is the only energy source for vegetation growth, and hence understanding the conversion efficiency of irradiance into biomass (light use efficiency) becomes essential for the simulation of NPP. Light use efficiency (LUE) is often assumed to be constant in the previous models (Whitehead et al., 2001). However, the daily variations in LUE did occur in New Zealand forest due to environmental limitations (Dugan & Whitehead, 2006). This becomes an indication for the improvement of NPP simulation model in this project.

This project uses Vensim software to design a simulation model of annual net primary production in a Wineberry and Fuschia dominated forest in the Taramakau Valley, Westland, New Zealand, which employs a bottom-up approach on the basis of the interactions between the forest growth and its available resources (including radiation, precipitation, and temperature).

2. Review of terrestrial net primary production modeling

Generally, net primary production models can be divided into two groups: empirical modeling and theoretically modeling. The first group models are based on the empirically-derived correlations of net primary productivity with available resources, and the other group of models are based on biochemical-mechanistic processes (Mulligan & Wainwright, 2004). In this section, both types of these models have been critically reviewed, with emphasis on their representative characteristics and limitations. Finally the choice of modeling types is discussed critically.

The Miami model, which empirically derives the associations between NPP and the mean annual temperature and precipitation, distinguishes plant types in the derivation of the correlation functions (Leith, 1975). This model can make reasonable estimations only for current NPP rates and their distribution, but not for the future NPP. This is mainly due to two reasons. Firstly, the Miami model does not consider the change of vegetation density. Secondly, the model assumption that vegetation distribution responses immediately to the climate change is sometimes not reasonable (Mulligan & Wainwright, 2004).

King et al., (1997) and Post et al., (1997) subsequently modified the Miami model. They quantified the changes in ecosystem carbon density by using a model of combining the empirical with the β-factor which added a carbon dioxide fertilization term to the Miami model. This β-factor was derived from the Farquhar photosynthesis formulation. A parameter (k1) was introduced as a scaling factor, which accounted for the efficiency losses in the actual increase of vegetation carbon density as a result of carbon dioxide fertilization. Particularly, this model distinguished relatively big amount of vegetation types in the global scope. However, the sub-model of β-factor did not take into account of the C4 photosynthesis pathway (King et al., 1997; Post et al., 1997).

Polglase and Wang (1992) also developed a biochemical model based on the β-factor. They distinguished ten types of biomes in this model. Baseline productivity was generated for each biome, which was based on the averages of the characteristic climate and productivity. Then β-factor modified the baseline productivity values. However, the parameters used in this β-factor were derived by using the average temperature. As the authors themselves suggested, significant inaccuracy would occur when there was significant variation in seasonal temperature (Polglase & Wang, 1992).

Lenton (2000) designed a quasi-biochemical model, in which photosynthesis was a function of temperature and atmospheric carbon dioxide, and respiration was depended on temperature and the vegetation carbon density. However, a particular problem can be caused by its assumption that the photosynthetic area is constant even though vegetation density may change (Mulligan & Wainwright, 2004).

Svirezhev and von Bloh (1997) developed another type of quasi-biogeochemical model. In this model, the maximum net primary productivity was generated as a parameter, which was modified by factors indicating the temperature response, water stress, carbon dioxide fertilization, and intra-specific competition. However, compared with other biochemical models, the functions interpreting the carbon dioxide fertilization and competition by Svirezhev and von Bloh (1997) were simply integrated into this model, which made it difficult to be understood (Mulligan & Wainwright, 2004).

Leuning et al., (1995) presented a process-based model. The model is based on the leaf-level measurements of photosynthesis and respiration, which is then scaled up to the canopy level. Sub-models of radiative transfer, energy balance, evaporation and photosynthesis are integrated. Stomatal conductance and photosynthesis is modified by a site water-balance sub-model. This model successfully predicts the forest canopy phenology (Leuning, et al., 1995). However, this model requires species-specific data for parameterization, which is not feasible for the application at a regional or global scale.

TRIFFID model calculates NPP by simulating photosynthesis and respiration. In this model, photosynthesis is calculated using the full Farquhar formulation, and respiration is determined by the rate of Rubisco carboxylation, both of which were modified by a water stress factor and the canopy level factor using Beer’s law. TRIFFID model is considered to be a relatively comprehensive vegetation model, incorporating intensively detailed information of biochemistry, physiognomy and forest dynamics. However, its complication leads it to be somewhat opaque and mathematically intractable (Cox, 2001; Mulligan & Wainwright, 2004).

Huntingford et al., (2000) subsequently presented a simplified version of the TRIFFID model, which was designed to investigate the qualitative response of the global vegetation to climate change. However, in this model, the description of a basic respiration required very detailed biochemical expression, and a single generic vegetation type relied on a broad parameterization, which still led this model to be inconvenient (Huntingford et al., 2000).

BIOME3 is one of the most comprehensive and complicated models for the estimation of NPP. It is also impenetrable and computationally intensive. In this model, photosynthesis is calculated by two simultaneous equations. The first relies on the Farquhar model modified by light and Rubisco limitations. An optimized daytime productivity and Rubisco carboxylation rates are given. The second is based on a stomatal conductance model that incorporates the effects of water stress. Then the NPP and leaf area index are adjusted by calculating an ‘extensive range’ of values from experiments. This model distinguishes different vegetation types by their maximum transpiration rates. Nevertheless, this complicated model requires a relatively high cost to apply (Haxeltine and Prentice, 1996; Mulligan & Wainwright, 2004).

The majority of these above models aim to simulate NPP in the global scale. However, because the environmental processes are usually scale-dependant, it is still not clear how these models interpret the data changes with the change of scales (Zhang & Wainwright, 2004).

Discussion of choice between empirical modeling and theoretical modeling: For the empirical models, the repeatable accuracy of modeling results, the availability of model input variables (such as convenience and economic cost accessing the modeling input data), the sensitivity of model parameters that influences the inaccuracy of modeling results (the more sensitive for the model parameters, the higher inaccuracy of modeling results, because all the parameters are estimated according to the empirical data, if the parameter is very sensitive, a tiny inaccuracy in parameter estimation will lead to a big error in final results) and the economic trade-off between the field-specific parameter determined by field sampling and estimation only of parameters are the main consideration factors to build a model; For the theoretical modeling, not only the empirical modeling results are assessed according to the above factors, but also the underlying theories are made more understandable. However, for the modeling of NPP at macro-scale, this article does not choose the theoretical biochemistry modeling as the first choice because biochemical modeling can only give better understanding at molecular scale rather than macro-scale, but these biochemical parameters increase the sensitivity of model parameters which may increase the inaccuracy of final modeling results.

3. Model description and implementation

This section will describe the equations used and their corresponding implementations in the Vensim system. Generally, the selection of these sub-equations facilitates the availability of input data. The final output of ‘Net primary production’ is calculated by an integrated equation. Among the five input variables, I employ the ‘Lookup’ function for the daily information input, while the ‘IF THEN ELSE’ function is used to facilitate the input process of quarterly data.

This model is based on a daily time scale with a spatial scale of m2. In the initial settings of this model, the initial time, final time, time step, and time unit is ‘1’, ‘365’, ‘1’, and ‘day’, respectively.

3.1. Absorbed Photosynthetically Active Radiation (PAR)

Fig 1. The interactions in sub-model of ‘Absorbed PAR’. (See PDF Article).

The daily absorbed PAR is considered as a function of ‘Daily radiation’ above canopy and ‘Daily fraction of absorbed PAR’. In the equation editing box of both ‘Daily radiation’ and ‘Daily fraction of absorbed PAR’, a ‘Lookup’ function is used for the daily data input. The unit of ‘Daily radiation’ is MJ/m2, while the ‘Daily fraction of absorbed PAR’ is without unit.

Fig 2. The editing box of ‘Lookup’ function for Daily radiation.(See PDF Article).

In the [Graph Lookup] box of ‘Daily radiation’, the days (x) axis is scale from 1 to 365, and the (y) axis is scaled by entering a value of 35 for the [Y-max]. Then the daily radiation data can be input in the left box (Fig 2).

Fig 3. The editing box of ‘Lookup’ function for ‘Daily fraction of absorbed PAR’.(See PDF Article).

Similarly, in the [Graph Lookup] box of ‘Daily fraction of absorbed PAR’, I scale the days (x) axis from 1 to 365, and (y) axis from 0 to 1. Then the data of daily fraction of absorbed PAR can be input in this box (Fig 3).

Fig 4. Equation editing for ‘Day’ variable.(See PDF Article).

Then a shadow variable <Time> is added, as a determinant to ‘Day’ (Fig 1). In the equation editing box of ‘Day’ variable, an equation of ‘Day = MODULO(Time, 365)+1’ is typed (Fig 4).

The variable of ‘Radiation’ is determined by both ‘Daily radiation’ and ‘Day’. Then, the variable of ‘Radiation’ can be interpreted as ‘Radiation = Daily radiation(Day)’ in the equation editing box (Fig 5). Similarly, for the variable ‘Fraction of PAR’, an equation of ‘Fraction of PAR = Daily fraction of absorbed PAR(Day)’ is added.

Fig 5. Equation editing for ‘Radiation’ variable.(See PDF Article).

However, not all the radiation can be utilized by forest canopy. To reflect this fact, I use ‘Light utilization’ linked with ‘Light saturation point’. The variable of ‘light saturation point’ is set to be a constant value (Q) of 22.04 MJ/m2. Then in the variable of ‘Light utilization’, an equation of‘IF THEN ELSE(Radiation<Light saturation point, Radiation,Light saturation point)’ is added.

Finally, the variable of ‘Absorbed PAR’ is determined by both ‘Fraction of PAR’ and ‘Light utilization’. The PAR is estimated to be 45% of total solar radiation. In the equation editing box, the equation is ‘Absorbed PAR = Fraction of PAR*Light utilization*0.45’.

3.2. Water limitation factor

Water deficit is one of the major factors limiting forest productivity. Hence, there will be a limitation factor of available water h(X) (Nicholas et al., (2007)) in this model, where X is the proportion of cumulative precipitation (∑P) to cumulative evapotranspiration (∑E) over a certain period. Parameter s is the slope of the curve relating saturation water vapor pressure to temperature, and γ is the psychrometric constant. Parameter a is empirical coefficient. In this model, the site-specific values of parameter a, s, and γ is estimated to be 1.8, 1.12, and 0.066, respectively, which are different from the value published in Nicholas et al., (2007) article.

‘Water limitation factor’ is determined by the input variables of ‘Radiation (3 months)’ and ‘precipitation (3 months)’ with a unit of ‘MJ/m2’ and ‘mm’, respectively. A shadow variable <Time> is added to determine these two input variables (Fig 6), which require quarterly climate data. For example, in 2007, the quarterly radiation values are 1500.24MJ, 511.28MJ, 730.27MJ, and 1598.99MJ from January to December at this site, and the 3-month precipitations are 282mm, 503.6mm, 369mm, and 532.4mm from January to December. Data are derived from NIWA.

Fig 6. The interactions in sub-model of water limitation factor.(See PDF Article).

I employ the IF THEN ELSE function to input these quarterly data. In the variable of ‘Radiation (3 months)’, an equation of‘IF THEN ELSE(Time<91, 1500.24, IF THEN ELSE(Time<182, 511.28, IF THEN ELSE(Time<274, 730.27, 1598.99)))’ is added (Fig 7). Similarly, the ‘Precipitation (3 months)’ is interpreted as ‘IF THEN ELSE(Time<91, 282, IF THEN ELSE(Time<182, 503.6, IF THEN ELSE(Time<274, 369, 532.4)))’. These two input variables determine the ‘Water balance’

Fig 7. The equation editing box for ‘Radiation (3 months)’.(See PDF Article).

Nicholas et al., (2007) equations contain two variables, the ‘Water balance’ and ‘Water limitation factor’ (Fig 6). In the equation editing box of ‘Water balance’, the equation of ‘1.8*"Precipitation (3 months)"/(1.12*0.45*"Radiation (3months)"/(1.12+0.066))’ is typed. I also use IF THEN ELSE function to edit ‘Water limitation factor’: ‘IF THEN ELSE(Water balance<1, Water balance, 1)’.

3.3. Light use efficiency

Dungan and Whitehead (2006) reported that atmospheric transmissivity explained most of daily variability in light use efficiency (ranging from 75% to 90%) for Fuchsia and Wineberry species in 1999 in this proposed forest. Therefore, in this project, light use efficiency (LUE) will be mainly considered as a function of atmospheric transimissivity (ε):

LUE = y0 + a × exp(-0.5×(ln(ε/x0)/b)2) equation 3

Based on their experiment results, the estimated values of parameter y0, a, x0, and b, are 0.28, 0.795, 0.18, and 0.78, respectively.

The atmospheric transmissivity is expressed as the fraction of irradiance reaching the top of the canopy after absorption and scattering by clouds and atmospheric turbidity (H/H0) (Bristow & Campbell, 1984). Bristow and Campbell (1984) proposed a model of estimating the atmospheric transmissivity (ε) which considered ε as a function of the difference between daily maximum and minimum temperature (ΔT). In this model, the parameters of aB, bB, cB are empirical coefficients with an estimated site-specific value of 0.66, 0.23, and 0.8, respectively, which is different from Bristow and Campbell (1984) article.

Fig 8. The interactions in sub-model of ‘Light use efficiency’.(See PDF Article).

A method which is similar to the ‘Daily radiation’ is employed for the input of ‘Daily temperature variation’. A ‘Lookup’ function is used for the ‘Daily T variation’ input. The shadow variable <Time> and a MODULO function are employed for the interpretation of variable ‘Days’. Then the ‘Temperature variation’ is determined by the ‘Daily T variation’ and ‘Days’ (Fig 8).

For the variable of ‘Atmospheric transmissivity’, the Bristow-Campbell equation ‘0.66*(1-EXP(-0.23*(Temperature variation^0.8)))’ is input in the equation editing box.

Then the ‘Light use efficiency’ is determined by ‘Atmospheric transmissivity’. The Dugan – Whitehead (2006) equation is typed in the equation editing box, which is ‘0.28+0.795*EXP(-0.5*((LN(Atmospheric transmissivity/0.18)/0.78)^2))’.

3.4. Net primary production

Finally, the ‘Net primary production’ will be calculated by an integrated equation:

NPP = PAR * FPAR * LUE * h(X) equation 5

Fig 9. The interaction between ‘Net primary production’ and its determinants.(See PDF Article).

In the equation box of ‘Net primary production’, a level variable is selected with an initial value of zero (Fig 10). Then the ‘Net primary production’ = INTEG [‘Absorbed PAR*Light use efficiency*Water limitation factor’] is typed in the below box.

Fig 10. The equation editing box of NPP.(See PDF Article).

4. Model test and validation

The model equations can be tested using the ‘check model’ function in the Vensim system. In this section, the output of each sub-model will be validated independently by comparing with either field studies or other model results.

4.1. Absorbed PAR

In the sub-model of absorbed PAR, the total estimated PAR is 1953.351MJ/m2 in 2007 (Fig 11). The mean fraction of absorbed PAR is 86.12%. Dugan and Whitehead (2006) did a field study in 1999 at the modeling site, and reported that the total PAR and the fraction of absorbed PAR were 3463MJ/m2 and 93%, respectively.

Fig 11. The absorbed PAR in 2007 in the modeling forest.(See PDF Article).

In this model, the PAR is estimated base on the assumption that the PAR is a constant proportion (45%) of total solar radiation. The total solar radiation data are measured by the meteorological station of Reefton Ews, which is approximately 6km away from the modeling forest. Surprisingly, the total estimated PAR in 2007 (1953.351MJ/m2) is significantly lower than the field measurements at the modeling site in 1999 (3463MJ/m2). In addition to the temporal and spatial variability in PAR, the different measurement methods used by Dugan and Whitehead (2006) may be also an important reason to explain this difference. In their field study, the PAR is directly measured by a quantum sensor (Li190SA, LiCor, Lincoln, Nebraska), which is different from the measurement equipment for total solar radiation used in the meteorological station.

The mean fraction of absorbed PAR in 2007 (86.12%) measured by Earth Resources Observation System Data Center (EROS) is also lower than the field measurement made by Dugan and Whitehead (2006) in 1999 (93%). Apart from the temporal variation in this fraction value, another reason was because the methods used by EROS calculated the absorbed PAR only and remove the PAR reflected by forest canopy, while Dugan and Whitehead (2006) calculated the intercepted PAR, which included both of them.

4.2. Light use efficiency

Fig 12. The estimated light use efficiency in 2007 in the modeling site.(See PDF Article).

The mean estimated light use efficiency by this model is approximately 0.65g C/MJ in 2007 (Fig 12). This result was consistent with Dugan and Whitehead (2006), who estimated that light use efficiency was approximately 0.60g C/MJ in 1999 in the modeling forest.

4.3. Water limitation factor

Fig 13. The estimated water limitation factor in 2007 in modeling site.(See PDF Article).

As can be seen in Fig 13, the productivity in 2007 was reduced by the water limitation factor only in the first 3 months. According to NIWA Climate database, the quarterly precipitation in 2007 was 282mm, 503.6mm, 369mm, and 532.4mm from January to December, which means that the first 3 months was the driest season in 2007. Therefore, the estimated value of water limitation factor matches to the precipitation records.

4.4. Net primary production

Fig 14. The estimated net primary production in 2007 in the modeling site.(See PDF Article).

The final output of net primary production in 2007 is 938.335g C/m2 (Fig 14). Dugan et al., (2004) conducted a simulation study using a process-based, biochemical-mechanistic model in the same forest in 2003, which showed that the estimated net primary production in 2003 was closely consistent with the estimation made by this project.

In addition, the net primary production of the modeling forest can be also obtained from the Earth Resources Observation System Data Center (EROS). According to their estimation, the average NPP was 806.7g C/m2 in 2007, which was slightly lower but still close to the result of this project. Therefore, the final output of 938.335g C/m2 in this project should be a reasonable value.

Compared with other global models, the estimated value of this project is lower than the simple TRIFFID, DEMETER, and TRIFFID models, but higher than BIOME3 model. Especially, the cumulative absorbed PAR in this project is 58 W/m2. This value reaches the light saturation point in TRIFFID model, above which little NPP increases with the increase of PAR in TRIFFID model.

5. Parameterization

There are nine parameters in total in this model. For some sub-models which have their own patent without publication of their parameterization methods, I used the reference values which are closest to the scenario of the modeling forest.

5.1. Atmospheric transmissivity

There are three empirical coefficients (aB, bB, cB) used in this sub-model. For the parameter aB, I use a reference value (0.66) calculated at a highland site which has similar topography to the proposed forest. For the bB and cB, Guillerdo et al., (2004) gave the equations for parameterization, which were calculated on the basis of the ΔT and φ, the daily difference between maximum and minimum temperature, and latitude of the proposed forest (42° 46' S), respectively. The ΔT data can be obtained from the NIWA Climate Database. I also built a simple model to calculate the bB, and cB values using Vensim software. The average calculated value of parameter bB and cB is 0.227 and 0.790, respectively. However, in this parameterization process, I used the daily ΔT data in 2007 only, which might not be sufficient to precisely determine these empirical coefficients (ideally, the data should cover the past ten years). Therefore, these are only approximate values.

5.2. Light use efficiency

Dugan and Whitehead (2006) proposed an equation (equation 3) to quantify the atmospheric transmissivity effect on the light use efficiency (LUE), without publication of parameter values. Nevertheless, they reported their relationship for each species (Fuchsia and Wineberry) (Fig 16).

Fig 16. The relationship between daily values of LUE and atmospheric transmissivity for fuchsiaand wineberry.Sources from Dugan and Whitehead (2006).(See PDF Article).

In the parameterization process, I firstly converted these species-specific results into mixed-species value, according to each species’ proportion of basal area in the modeling forest. This conversion equation is:

LUE (mixed species) = LUE (Fuchsia)*0.063 + LUE (Wineberry)*0.937

Then four paired results, which are (0.1, 0.43), (0.2, 0.62), (0.4, 0.49), and (0.6, 0.38), are obtained for the calculation of these parameter values in equation 3. The estimated parameter values of y0, a, x0, and b are 0.28, 0.795, 0.18, and 0.78, respectively.

5.3. Light saturation point

The light saturation point (Q) used in this model was estimated from the global model TRIFFID (Fig 13). This model result showed a limitation value yearly for light utilization.However, I converted this yearly value into a daily value of 20.02MJ/m2.

However, this light saturation point which is based on a daily temporal scale sometimes may not properly explain the photosynthesis process (e.g. light saturation may occur when the daily radiation is below this limitation value if the hourly radiation is high during that day). Most physiological studies employ a time scale of seconds to determine the light saturation value. In practice, the NPP shows an asymptotic relationship with the absorbed PAR. However, this parameterization of Q assumes that the relation between NPP and PAR is linear.

5.4. Water limitation factor

In equation 1, parameter ‘s’ is the slope of the curve relating saturation water vapor pressure to temperature, and ‘γ’ is the psychrometric constant. These two parameters are constant, which can be found in the references (“HyperPhysics: Thermodynamics, and The Asce Standardized Reference Evapotranspiration Equation”). However, the empirical coefficient ‘a’ is site-specific. I used a reference value of 1.8 which was derived from a temperate Douglas-fir forest in Australia.

6. Sensitivity analysis

Sensitivity analysis not only contributes to the model validation, but also becomes essential for the future model development. There are a number of reasons to do sensitivity analysis for the modeling. Firstly, sensitivity analysis helps to better understand the behavior and interactions between the output and parameters; in the multi-component model, the sub-models can verified by sensitivity analysis; it helps to improve the model parsimony by rejecting the insensitive parameters; it contributes to identify which parameters require additional field research to strengthen the knowledge base; sensitivity analysis also examines which inputs are highly correlated to the final outputs (Mulligan & Wainwright, 2004).

Hamby, (1994) summarized a large number of sensitivity analysis methods for model parameters. In this project, the sensitivity analysis was implemented using Pearson’r correlation methods. Ten tests in total were employed for this correlation test. In each test, all the parameters were randomly selected within a reasonable range which was determined by historical records. Then the correlation (R) between each parameter and final output was calculated using IFA Online Statistic Service (Table 1). The higher the correlation (R), the more sensitive is the parameter.

Table 1. Sensitivity analysis for nine parameters in this model.(See PDF Article).

However, there are some drawbacks in this correlation test used for sensitivity analysis: the method is based on the assumption that the input/output relationship is linear (β is an exponential parameter in this model); and some input parameters strongly correlated to one another (such as aB and β) may result in significant input/output correlations (Hamby, 1994).

In each test, parameter aB, bB, cB, y0, β, x0, b, a, and Q are random selected within a range of [0.5, 0.9], [0.15, 0.30], [0.70, 1.00], [0.15, 0.40], [0.40, 0.70], [0.10, 0.30], [0.70, 1.00], [1.0, 4], and [18.00, 25.00], respectively. As can be seen from the results (Table 1), parameter aB, β, x0, and a are highly correlated with the final output, which show a high sensitivity in this model.

Particularly, the parameter Q is insensitive in this model. As discussed above, this light saturation point (Q) does have limitation in this model. Therefore, this parameter Q may be rejected to ensure the parsimony of the model, or requires further studies in the future development.

7. Conclusion

7.1. Application and Advantages

This simulation model well estimated the net primary production in 2007, which has been validated in section 4.4. Among these five input variables, daily precipitation, daily radiation, 3-month radiation, 3-month precipitation, and daily variation of temperature can be downloaded from NIWA Climate Database. The fraction of absorbed PAR can be derived from Earth Resources Observation System Data Center (EROS), with a spatial scale of 1 km2. These fraction data are processed per 8 days, which can be downloaded in HDF format, and then can be transferred into excel for the calculation of mean values. Therefore, all of the data used in this model is readily and freely available for forest managers. Another advantage of this model is that it well indicates the daily variations in LUE, which improves the accuracy of simulating forest growth.

Another advantage of this numerical modeling is that the parameters can be estimated on the basis of empirical values to reduce the modeling cost. For example, the parameters of water limitation factor or light saturation point can be estimated according to the past study on the relevant species, which have been effectively done in my study. There are four steps in parameterization: estimation of parameters according to the empirical values; conducting sensitivity analysis for evaluating the accuracy and credibility of estimated parameters; validation of the estimated parameters by comparing the final modeling results with the field observations at experimental scale; final promotion of this modeling in larger scale. Generally, the more site-specific or species-specific parameters required by field measurement, the higher cost for building a model. Then there is a trade-off between the field measurement and the estimation based on the empirical value, which have been deduced in my modeling study.

7.2. Limitations and future development

7.2.1. Lack of considerations on the available soil nutrient and pest impact

In practice, soil available nutrients, such as nitrogen, potassium, and phosphate, usually become a limitation for the forest growth. However, this model is based on the assumption that there are always sufficient soil nutrients in this proposed forest. In addition, this model also assumes that the proposed forest is pest-free, which contradicts the fact that New Zealand forest usually suffers from a number of disease, insects, and invasive mammals.

7.2.2. Temporal and spatial scale

This model employed a daily time scale, which matches to the meteorological data from NIWA Climate Database. It may not be feasible to predict the future annual NPP using this model, because it seems to be impossible to predict the future long-term climate data on the basis of a daily scale.

Lots of environmental process and patterns are scale-dependent. In this model, the m2 of spatial scale is used, which is also a common spatial scale in other NPP estimation models. The empirical data used for parameterization in this model is based on the local scale (the whole forest range), and then converted into 1 m2. Nevertheless, the relatively small scale makes it impossible to indicate the species composition and distribution. Further more, this model cannot be applied to other forests in NZ with different species composition, which requires additional field studies for the empirical parameterization.

7.2.3. Future model development

This model is non-spatial, which does not consider the ‘species flow’. However, the ‘Light use efficiency’ used in this model can be significantly influenced by the change of species composition, because variations in LUE do occur among different species. Further more, the application of 3S technology on the biodiversity monitoring of vegetation ecosystem, novelly presented by the article 3 of this journal, facilitates the integration of species biodiversity on this carbon sink model in global scale.

In the parameterization process, I integrate the species-specific LUE of two species into a ‘mixed-species’ value according to their proportion of basal areas in the proposed forest. However, this proportion may change with the change of species composition over a longer period. The future model development should integrate the species composition change into the determination of ‘mixed-species’ LUE.

This is the revised materials in book “Proceedings for Degree of Postgraduate Diploma in Environmental Science (3rd Edition).” Published in 2016. The ‘chapter’ content mentioned in this article is in previous book. Revised on 25/01/2021. Secondly Revised on 01/09/2021. Thirdly Revised on 09/01/2022. Fourthly Revised on 29/11/2022. Fifthly revised on 24/04/2023. Latest revised on 31/05/2023.

Adams, B., White, A., & Lenton, T. M. (2004). An analysis of some diverse approaches to modelling terrestrial net primary productivity. Ecological Modelling [Ecol. Model.]. 177, 3-4. https://doi.org/10.1016/j.ecolmodel.2004.03.014

Akihiko, I. A., & Takehisa, O. (2002). A simulation model of the carbon cycle in land ecosystems (Sim-CYCLE): a description based on dry-matter production theory and plot-scale validation. Ecological Modelling [Ecol. Model.]. 151, 143-176. https://doi.org/10.1016/S0304-3800(01)00473-2

Cleugh, H. A., Leuning, R., Mu, Q., & Running, S. W. (2007). Regional evaporation estimates from flux tower and MODIS satellite data. Remote Sensing of Environment [Remote Sens. Environ.]. 106(3), 285-304. https://doi.org/10.1016/j.rse.2006.07.007

Conover, W. J. (1980). Practical Nonparametric Statistics. New York: John Wiley & Sons.

Dungan, R. J., & Whitehead, D. (2006). Modelling environmental limits to light use efficiency for a canopy of two broad-leaved tree species with contrasting leaf habit. New Zealand Journal of Ecology [N. Z. J. Ecol.]. 30(2), 251-259.

Dungan, R. J., Whitehead, D., Mcglong, R., Allen, B., & Ducan, R. (2004). Simulated carbon uptake for a canopy of two broadleaved tree species with contrasting leaf habit. Functional Ecology [Funct. Ecol.]. 18, 34-42. https://doi.org/10.1111/j.1365-2435.2004.00809.x

Garcia, J. V. (1994). Principios F'ısicos de la Climatolog'ıa.: Ediciones UNALM.

Guillermo, A., Baiggrio, A., Esequiel, B., Villegas, C., Irene, T., Jose, F., et al. (2004). Atmospheric transmissivity: distribution and empirical estimation around the central andes. International Journal of Climatology, 24, 1124-1136. https://doi.org/10.1002/joc.1060

Hall, F. G., Huemmerich, K. F., & Goward, S. N. (1990). Use of narrow-band spectra to estimate the fraction of absorbed photosynthetically radiation. Remote Sensing and Environment, 15, 271-282. https://doi.org/10.1016/0034-4257(90)90097-6

Hamby, D. M. (1994). A review of methodology for parameter sensitivity analysis in environmental modelling. Environmental monitoring and assessment., 32, 135-154. https://doi.org/10.1007/BF00547132

Hargreaves, G., & Samani, Z. (2002). Estimating potential evapotranspiration. Journal of Irrigation and Drainage Engineering, 108, 225-230. https://doi.org/10.1061/JRCEA4.0001390

Helton, J. C., Garner, J. W., Marietta, M. G., Rechard, R. E., Rudeen, D. K., & Swift, E. (1993). Uncertainty and Sensitivity Analysis Results Obtained in a Preliminary Performance assessmentfor the Waste Isolation Pilot Plant. Nuc. Sci. and Eng, 114, 286-331. https://doi.org/10.13182/NSE93-A24041

Hoffman, E., & Gardner, R. H. (1983). 'Evaluation of Uncertainties in Environmental RadiologicalAssessment Models'. In H. R. Meyer (Ed.), Radiological Assessments: a Textbook on Environmental Dose Assessment. Washington, DC: Nuclear Regulatory Commission Report No. NUREG/CR-3332.

Keith L. Bristow, Gaylon S. Campbell (1984),On the relationship between incoming solar radiation and daily maximum and minimum temperature, Agricultural and Forest Meteorology, Volume 31, Issue 2, 1984, Pages 159-166, ISSN 0168-1923, https://doi.org/10.1016/0168-1923(84)90017-0. https://doi.org/10.1016/0168-1923(84)90017-0

Leuning, R., Helen, A., Cleugh, S. J., & Zegelin, D. H. (2005). Carbon and water fluxes over a temperate Eucalyptus forest and a tropical wet/dry savanna in Australia: measurements and comparison with MODIS remote sensing estimates. Agricultural and Forest Meteorology [Agric. For. Meteorol.]. 129, 151-173. https://doi.org/10.1016/j.agrformet.2004.12.004

Manzul, K. H., & Ysuoka, Y. (2005). Estimation of net primary productivity by integrating remote sensing data with an ecosystem model. Remote Sensing of Environment [Remote Sens. Environ.]. 94, 298-310. https://doi.org/10.1016/j.rse.2004.10.004

Mu, Q., Heinsch, F. A., Running, S. W., Cleugh, H. A., & Leuning, R. (2006). Evaluation of Satellite-based Evapotranspiration Using Fluxnet Data. 2000 Florida Ave., N.W. Washington DC: American Geophysical Union.

Mulligan, M., & Wainwright, J. (2004). Modelling and Modelling Building. In J. Wainwright & M. Mulligan (Eds.), Environmental Modelling: Finding Simplicity in Complexity: Wiley.

Myneni, P., Hoffman, B., Knyazikhin, A., Privette J L, J, G., Tian Y, et al. (2002). Global products of vegetation leaf area and fraction absorbed PAR from year one of MODIS data. Remote Sensing of Environment [Remote Sens. Environ.]. 83, 214-231. https://doi.org/10.1016/S0034-4257(02)00074-3

Nicholas, C., Coops, A., Rachhpal, S., Jassal, B., Leuning R C, Andy, T., et al. (2007). Incorporation of a soil water modifier into MODIS predictions of temperate Douglas-fir gross primary productivity: Initial model development. Agricultural and Forest Meteorology [Agric. For. Meteorol.]. 147, 99-109. https://doi.org/10.1016/j.agrformet.2007.07.001

Polglase, P. J., & Wang, Y. P. (1992). Potential CO2-enhanced Carbon Storage by the Terrestrial Biosphere. Aust. J. Bot., 40, 641-656. https://doi.org/10.1071/BT9920641

Post, W. M., King, A. W., & Wullshleger, S. D. (1997). Historical variations in terrestrial biospheric carbon storage. Global Biogeochem. Cycles, 11(1), 99-109. https://doi.org/10.1029/96GB03942

Rosati, A., Metcalf, S. D., & Lampine, B. D. (2004). A Simple Method to Estimate Photosynthetic Radiation Use Efficiency of Canopies. Annals of Botany, 93, 567-574. https://doi.org/10.1093/aob/mch081

Smith, B., Knorr, W., Widlowski, J. L., Pinty, B., & Gobron, N. (2008). Combining remote sensing data with process modelling to monitor boreal conifer forest carbon balances. Forest Ecology and Management [For. Ecol. Manage.]. 255(12), 3985-3994. https://doi.org/10.1016/j.foreco.2008.03.056

Turner, D. P., Urbanski, S., Bremer, D., Wofsy, S. C., Meyers, T., Gower, S. T., et al. (2003). A cross-biome comparison of daily light use efficiency for gross primary production. Global Change Biology, 9, 383-395. https://doi.org/10.1046/j.1365-2486.2003.00573.x

Zimmermann, H. J. (2000). An application-oriented view of modelling uncertainty. European Journal of Operational Research, 122, 190-198. https://doi.org/10.1016/S0377-2217(99)00228-3