CN108764583B - Unbiased prediction method for forest accumulation - Google Patents

Unbiased prediction method for forest accumulation Download PDF

Info

Publication number
CN108764583B
CN108764583B CN201810575883.7A CN201810575883A CN108764583B CN 108764583 B CN108764583 B CN 108764583B CN 201810575883 A CN201810575883 A CN 201810575883A CN 108764583 B CN108764583 B CN 108764583B
Authority
CN
China
Prior art keywords
model
forest
parameter
noise
order
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810575883.7A
Other languages
Chinese (zh)
Other versions
CN108764583A (en
Inventor
刘恩斌
周国模
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang A&F University ZAFU
Original Assignee
Zhejiang A&F University ZAFU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang A&F University ZAFU filed Critical Zhejiang A&F University ZAFU
Priority to CN201810575883.7A priority Critical patent/CN108764583B/en
Publication of CN108764583A publication Critical patent/CN108764583A/en
Application granted granted Critical
Publication of CN108764583B publication Critical patent/CN108764583B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/02—Agriculture; Fishing; Forestry; Mining

Landscapes

  • Business, Economics & Management (AREA)
  • Engineering & Computer Science (AREA)
  • Strategic Management (AREA)
  • Human Resources & Organizations (AREA)
  • Economics (AREA)
  • Theoretical Computer Science (AREA)
  • Marketing (AREA)
  • General Physics & Mathematics (AREA)
  • General Business, Economics & Management (AREA)
  • Physics & Mathematics (AREA)
  • Tourism & Hospitality (AREA)
  • Animal Husbandry (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Agronomy & Crop Science (AREA)
  • Mining & Mineral Resources (AREA)
  • General Health & Medical Sciences (AREA)
  • Marine Sciences & Fisheries (AREA)
  • Health & Medical Sciences (AREA)
  • Primary Health Care (AREA)
  • Development Economics (AREA)
  • Game Theory and Decision Science (AREA)
  • Entrepreneurship & Innovation (AREA)
  • Operations Research (AREA)
  • Quality & Reliability (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

The invention discloses a forest accumulation amount estimation and prediction method, which is sequentially carried out according to the following steps: firstly, collecting time series paired data of an object region; secondly, determining the order of the autoregressive moving average model; thirdly, solving a primary estimator of the parameter; fourthly, solving a noise model parameter vector; fifthly, calculating a filtering signal; sixthly, solving secondary estimation quantity of the parameters; and seventhly, calculating the final parameter estimation value to obtain an estimation value and a predicted value of the forest accumulation. The method for estimating and predicting the forest accumulation has the advantages of five aspects: factors which can not be obtained by comprehensive actual measurement and influence the accumulation amount can be comprehensively measured; the method can describe the autocorrelation of the accumulation amount and reflect the regression relationship between the accumulation amount and the forest area; the nonlinearity of the tree growth can be reflected; the method has a complete theoretical basis and can overcome the defects of an empirical model; by introducing the whitening filter, unbiased estimation of the model parameters can be obtained, and the defect of biased estimation of the common least square can be overcome.

Description

Unbiased prediction method for forest accumulation
Technical Field
The invention relates to a forest accumulation amount estimation and prediction method, in particular to a forest amount estimation and prediction method based on an autoregressive moving average model and generalized least square.
Background
The forest accumulation is an important component of regional natural resources and is one of the most important indexes for forest resource investigation, monitoring and prediction. The method has the advantages that the forest accumulation is estimated and predicted efficiently, quickly and accurately, the basis for estimating the carbon budget of the forest is realized, and a decision basis can be provided for sustainable development of forest resources and strategic formulation of regional forest development. The forest accumulation amount is quickly and accurately estimated and predicted, the method belongs to the scientific category oriented to decision making, and the estimation and prediction method is important.
The current methods for estimating and predicting forest accumulation mainly comprise four methods: the method has the advantages that the method is suitable for simulating a complex system, has the functions of self-learning, associative storage and high-speed searching of optimal solutions, and has the main defects that the convergence speed is low, the number of hidden layers and hidden nodes is difficult to determine, and the estimation accuracy is influenced by the convergence speed; secondly, a grey system theory, namely, the method has the advantages that firstly, the satisfactory prediction effect is difficult to obtain on the original non-stationary time sequence, and secondly, the known conditions lack strict theoretical basis, so that the estimation result is not necessarily the best prediction; thirdly, a K-nearest neighbor method (K-NN) is used for researching forest accumulation of Italy two regions, and the research shows that estimation and prediction results are influenced by factors such as remote sensing data spectrum types, input auxiliary variables, types of measured multi-dimensional distances, number of nearest neighbor samples and the like; fourthly, the application of 3S technology, the method has the advantages that a fast, economic, convenient and reliable way is provided for large-scale forest accumulation estimation and prediction, the defect is that remote sensing data has limitations, especially in tropical and subtropical areas with cloud rain and fog, visible light and infrared remote sensing are greatly limited, factors such as forest age, canopy density, soil thickness, dominant tree species and the like are generally selected for accumulation estimation and prediction variables, but the factors are difficult to obtain through remote sensing and GIS, in addition, polynomial models are mostly selected for accumulation estimation and prediction models, as forest distribution, types, growth conditions and the like are large, and the polynomial equation is an empirical model and has great subjectivity, the polynomial models are not necessarily suitable for forest accumulation estimation and prediction.
The existing forest accumulation amount prediction and estimation method is greatly influenced by noise pollution in the calculation process, and the result always has small deviation. Therefore, the problem that technical personnel in the industry continue to solve is to design a forest accumulation estimation and prediction method which can weaken or even eliminate noise pollution, has accurate result, can fully utilize information of forest resource survey data, objectively reflect the total scale and the abundance of forest resources, and further can measure the quality of forest ecological environment and provide reliable data support for regional forest development planning.
Disclosure of Invention
The invention aims to provide an unbiased estimation method for forest accumulation. The method effectively weakens or even eliminates the influence of noise pollution, has more accurate result, can fully utilize the information of forest resource survey data, and can objectively reflect the total scale and the abundance degree of the forest resource, thereby measuring the quality of the forest ecological environment and providing reliable data support for regional forest development planning.
The technical scheme of the invention is as follows: the unbiased prediction method of the forest accumulation amount comprises the following steps:
1) collecting the paired data of the latest continuous p-period forest land area and forest accumulation amount actual measurement time sequence of the target region, wherein p is more than 2n +1, and n is the order of the model in the formula (1);
2) estimating the forest accumulation amount of the target region based on the autoregressive moving average model and the generalized least square algorithm through the data collected in the step 1), wherein the method specifically comprises the following steps:
a. a forest accumulation time sequence model and a noise lag operator polynomial:
the time series model of the forest accumulation is set as follows:
Figure BDA0001686893020000021
wherein: n is the order of the model; k. j is constant and k is more than j, y (k) is the output of the system and represents the forest accumulation amount of the kth period; y (k-j) is the historical output of the system, representingForest accumulation amount in the k-j period; u (k-j) is a system determinability input and represents the forest land area at the k-j stage; a isjIs an autoregressive coefficient, bjFor input of transfer coefficient, i.e. ajAnd bjParameters of a forest accumulation time series model; ξ (k) is the residual error of the kth stage, i.e. the colored noise of the model;
let the noise lag operator polynomial be:
Figure BDA0001686893020000031
wherein: z is a radical of-i(i ═ 1, 2.. times.n) is a hysteresis operator, a1,a2,…,an,b0,b1,b2,…,bnIs a parameter;
equation (1) can be expressed as:
a(z-1)y(k)＝b(z-1)u(k)+ξ(k) (3)
b. determining the order n of the model:
let θ be (a)1,a2,…,an,b0,b1...,bn) Representing model parameter vectors;
performing least square fitting by using models of different orders to calculate goodness of fit
Figure BDA0001686893020000032
Wherein
Figure BDA0001686893020000033
Is the parameter θ of the model at a certain orderTEstimating quantity;
when the order of the model is increased, the model is,
Figure BDA0001686893020000034
the value of (a) is reduced, and when the significant reduction is terminated, the order of the corresponding model is the appropriate model order;
c. calculating a first estimate of the parameter, i.e. the quantity that the generalized least squares algorithm needs to solve first:
let the whitening filter, i.e. the whitening function, be:
f(z-1)＝1+f1z-1+...+fmz-m； (4)
wherein: z is a radical of-j(j ═ 1, 2.. times, m) is a hysteresis operator, f1,f2,…,fmIs a parameter;
the colored noise ξ (k) is transformed into white noise w (k) as follows:
Figure BDA0001686893020000035
wherein: f (z)-1) Is a whitening filter;
let f (z) in formula (5)-1) 1, i.e. the noise model parameter vector f is an m-dimensional zero vector, where f ═ f1,f2,...,fm) Then model a (z)-1)y(k)＝b(z-1) u (k) + ξ (k) or
Figure BDA0001686893020000041
The residual error of (a) is degraded to white noise, i.e., ξ (k) ═ w (k);
order to
Figure BDA0001686893020000042
Wherein:
Figure BDA0001686893020000043
information is output for the n-period history,
Figure BDA0001686893020000044
inputting information for the n +1 period;
when p pairs of sampling values { x ] are obtainedT(k) Y (k) }, k ═ n +1, n + 2.., n + p, according to y (k) ═ xT(k)θT+ ξ (k), the vector matrix model is derived:
Y＝XθT+ξ； (6)
wherein: y isT＝(y(n+1),y(n+2),...,y(n+p))，ξT＝(ξ(n+1),ξ(n+2),...,ξ(n+p))，
Figure BDA0001686893020000045
Determining theta by ordinary least square methodTFirst approximation estimator of (1):
Figure BDA0001686893020000046
d. and (3) solving a noise model parameter vector:
order to
Figure BDA0001686893020000047
φT(k) Representing a colored noise sequence under random interference of a forest ecosystem; by using
Figure BDA0001686893020000048
And paired samples of each period, consisting of y (k) xT(k)θT+ ξ (k) the estimated value of colored noise ξ (k) at each stage of the model is calculated:
Figure BDA0001686893020000049
will obtain
Figure BDA00016868930200000413
Respectively substitute for
Figure BDA00016868930200000410
Obtaining:
Figure BDA00016868930200000411
according to the vector matrix model:
Figure BDA00016868930200000412
wherein:
Figure BDA0001686893020000051
Figure BDA0001686893020000052
using the common least square method to obtain fTThe estimated amount of (c):
Figure BDA0001686893020000053
wherein:
Figure BDA0001686893020000054
is composed of
Figure BDA0001686893020000055
Transposing;
e. calculating a filtered signal:
substituting formula (5) into a (z)-1)y(k)＝b(z-1) u (k) + xi (k) to obtain
a(z-1)f(z-1)y(k)＝b(z-1)f(z-1)u(k)+w(k)；
Using the products obtained in step d
Figure BDA0001686893020000056
Calculating the input filtered signal:
Figure BDA0001686893020000057
and outputting a filtered signal:
Figure BDA0001686893020000058
and input and output filtered signals:
Figure BDA0001686893020000059
f. solving a secondary estimator of the parameters:
according to step d, finding
Figure BDA00016868930200000510
Then according to step e, obtain
Figure BDA00016868930200000511
And
Figure BDA00016868930200000512
by using
Figure BDA00016868930200000513
And
Figure BDA00016868930200000514
respectively in place of the formula f (z)-1)y(k)＝f(z-1)xT(k)θTY (k) and x in + w (k)T(k) And obtaining a common linear regression model between the output filtering signal and the input and output filtering signal:
Figure BDA00016868930200000515
when k is n +1, n +2, …, n + p, the vector matrix model is obtained:
Figure BDA00016868930200000516
wherein:
Figure BDA00016868930200000517
Figure BDA0001686893020000061
determining theta by using ordinary least square methodTSecond order approximation estimator of (1):
Figure BDA0001686893020000062
g. calculating the final parameter estimator:
by using
Figure BDA0001686893020000063
Substitution
Figure BDA0001686893020000064
According to step d, calculating
Figure BDA0001686893020000065
According to step f, calculating to obtain a third approximate estimator of the parameter
Figure BDA0001686893020000066
Repeating the iteration until
Figure BDA0001686893020000067
Wherein: i is
Figure BDA0001686893020000068
N is
Figure BDA0001686893020000069
The total number of iterations of (c);
at this time, θTIs estimated by
Figure BDA00016868930200000610
Then will be
Figure BDA00016868930200000611
Substitution into
Figure BDA00016868930200000612
And obtaining an estimated value or a predicted value of y (k), namely obtaining an estimated value or a predicted value of the forest accumulation amount of the k-th period of the target region.
Compared with the prior art, the method and the device can accurately estimate and predict the forest accumulation amount by starting from the following two aspects:
one is time series analysis. The method is a widely applied data analysis method and is mainly used for describing and exploring the quantitative rule of the phenomenon changing along with the development of time. An important aspect of time series analysis is to estimate the current situation and predict the future based on the data already in the past. Because factors influencing forest accumulation are complicated and data information of some influencing factors cannot be obtained, the time series analysis method for comprehensively substituting the factors by time t can achieve the purposes of estimation and prediction. The autoregressive moving average model, namely the ARMA model, is established according to the self measuring factor time sequence of forest growth without considering the operation condition in an ecological system, so that the ARMA model is a better model for solving a nondeterministic system, can not only make state estimation on the system, but also can perform trend prediction on the operation condition in the system. Further, the forest accumulation is a variable related to time and area of the wooded land, so the preferred model for estimating and predicting forest accumulation is an ARMA model of time series analysis. As the forest land area and the forest accumulation amount of each period of resource investigation are paired sampling data, an ARMA (n, n) model is selected;
and secondly, identifying each parameter of the ARMA model by using a generalized least square algorithm. The most common parameter identification method of the ARMA model is a common least square algorithm, however, the common least square algorithm is unbiased only when the residual term of the ARMA model is white noise. Due to random interference of a plurality of natural and unnatural factors, the internal characteristics of the forest ecosystem have randomness, so that the input end of the ARMA model is inevitably polluted by noise, and the residual term of the ARMA model is colored noise in general. Therefore, the estimation and prediction of forest accumulation by using the ordinary least square algorithm is biased. The generalized least square algorithm filters colored noise into white noise, namely, ξ (k) is converted into w (k), so that an unbiased estimation value of the forest accumulation is obtained, and the estimation and prediction precision of the forest accumulation is improved.
The model of this application can fully reflect the growth law of forest self and relatively accurate. The forest ecosystem is a random system, the accumulation amount of the forest ecosystem is related to a plurality of factors, accumulation or consumption of accumulation in the next period is accumulation or consumption in the first period, so that a certain autocorrelation relation exists between the accumulation, influence accumulation factors which cannot be obtained through actual measurement are integrated by using time t, and an ARMA accumulation model related to the forest land area can reflect the change rule of the accumulation. The common forest accumulation estimation method generally adopts a method of a deterministic and constant system to predict forest dynamics, which is not consistent with the characteristics of forest resource change, and the influence of random interference on a prediction model is considered, and the forest change is described by using a random variable, so that the result is more accurate.
The forest growth sequence is a non-stationary sequence and therefore it is relatively complex and practically difficult to describe this process accurately. The autoregressive moving average model is a dynamic model taking the state characteristics of the forest as the dominant information, and can accurately describe the accumulation or consumption process of forest growth in a non-stationary state. Meanwhile, the model has a complete theoretical basis and can overcome the defect that some estimation methods adopt empirical models. The autoregressive moving average model can also be regarded as a systematic difference equation model, and the difference equation can describe some nonlinear models conforming to the growth of trees, such as let y (t) etThen, then
Figure BDA0001686893020000084
Let y (t) ta(a is an integer), then
Figure BDA0001686893020000081
Wherein
Figure BDA0001686893020000082
Is a parameter, thereby (t +1)aCan also be expressed by a difference equation. The model is combined with generalized least square to obtain unbiased estimation of parameters, and the defect of biased estimation of common least square can be overcome.
In conclusion, the method effectively weakens or even eliminates the influence of noise pollution, has more accurate result, can fully utilize the information of forest resource survey data, and can objectively reflect the total scale and the abundance of the forest resources, thereby measuring the quality of the forest ecological environment and providing reliable data support for regional forestry development planning.
Detailed Description
The present invention is further illustrated by the following examples, which are not to be construed as limiting the invention.
Example (b): the unbiased prediction method of the forest accumulation amount comprises the following steps:
1) collecting the paired data of the latest continuous p-period forest land area and forest accumulation amount actual measurement time sequence of the target region, wherein p is more than 2n +1, and n is the order of the model in the formula (1). When paired data is collected, the more the period number p is, the better, and the value of n is within an estimated range (such as n is 2 or n is 3), so that the collected period number p basically satisfies p > 2n + 1.
2) Estimating the forest accumulation amount of the target region based on the autoregressive moving average model and the generalized least square algorithm through the data collected in the step 1), wherein the method specifically comprises the following steps:
a. a forest accumulation time sequence model and a noise lag operator polynomial:
the time series model of the forest accumulation is set as follows:
Figure BDA0001686893020000083
wherein: n is the order of the model; k. j is constant and k is more than j, y (k) is the output of the system and represents the forest accumulation amount of the kth period; y (k-j) is the historical output of the system and represents the forest accumulation amount in the k-j stage; u (k-j) is a system determinability input and represents the forest land area at the k-j stage; a isjIs an autoregressive coefficient, bjIn order to input the transfer coefficient, the transfer coefficient is input,namely ajAnd bjParameters of a forest accumulation time series model; ξ (k) is the residual error of the kth stage, i.e. the colored noise of the model;
let the noise lag operator polynomial be:
Figure BDA0001686893020000091
wherein: z is a radical of-i(i ═ 1, 2.. times.n) is a hysteresis operator, a1,a2,…,an,b0,b1,b2,…,bnIs a parameter;
equation (1) can be expressed as:
a(z-1)y(k)＝b(z-1)u(k)+ξ(k) (3)
b. determining the order n of the model:
let θ be (a)1,a2,...,an,b0,b1...,bn) Representing model parameter vectors;
performing least square fitting by using models of different orders to calculate goodness of fit
Figure BDA0001686893020000092
Wherein:
Figure BDA0001686893020000093
is the parameter θ of the model at a certain orderTEstimating quantity;
when the order of the model is increased, the model is,
Figure BDA0001686893020000094
the value of (a) is reduced, and when the significant reduction is terminated, the order of the corresponding model is the appropriate model order;
c. calculating a first estimate of the parameter, i.e. the quantity that the generalized least squares algorithm needs to solve first:
let the whitening filter, i.e. the whitening function, be:
f(z-1)＝1+f1z-1+…+fmz-m； (4)
wherein: z is a radical of-j(j ═ 1, 2.. times, m) is a hysteresis operator, f1,f2,…,fmIs a parameter;
the colored noise ξ (k) is transformed into white noise w (k) as follows:
Figure BDA0001686893020000095
wherein: f (z)-1) Is a whitening filter;
let f (z) in formula (5)-1) 1, i.e. the noise model parameter vector f is an m-dimensional zero vector, where f ═ f1,f2,...,fm) Then model a (z)-1)y(k)＝b(z-1) u (k) + ξ (k) or
Figure BDA0001686893020000101
The residual error of (a) is degraded to white noise, i.e., ξ (k) ═ w (k);
order to
Figure BDA0001686893020000102
Wherein:
Figure BDA0001686893020000103
information is output for the n-period history,
Figure BDA00016868930200001015
inputting information for the n +1 period;
when p pairs of sampling values are obtained
Figure BDA0001686893020000104
When, according to y (k) ═ xT(k)θT+ ξ (k), the vector matrix model is derived:
Y＝XθT+ξ； (6)
wherein: y isT＝(y(n+1),y(n+2),...,y(n+p))，ξT＝(ξ(n+1),ξ(n+2),...,ξ(n+p))，
Figure BDA0001686893020000105
Determining theta by ordinary least square methodTFirst approximation estimator of (1):
Figure BDA0001686893020000106
d. and (3) solving a noise model parameter vector:
order to
Figure BDA0001686893020000107
φT(k) Representing a colored noise sequence under random interference of a forest ecosystem; by using
Figure BDA0001686893020000108
And paired samples of each period, consisting of y (k) xT(k)θT+ ξ (k) the estimated value of colored noise ξ (k) at each stage of the model is calculated:
Figure BDA0001686893020000109
will obtain
Figure BDA00016868930200001010
Respectively substitute for
Figure BDA00016868930200001011
Obtaining:
Figure BDA00016868930200001012
according to the vector matrix model:
Figure BDA00016868930200001013
wherein:
Figure BDA00016868930200001014
Figure BDA0001686893020000111
using the common least square method to obtain fTThe estimated amount of (c):
Figure BDA0001686893020000112
wherein:
Figure BDA0001686893020000113
is composed of
Figure BDA0001686893020000114
Transposing;
e. calculating a filtered signal:
substituting formula (5) into a (z)-1)y(k)＝b(z-1) u (k) + xi (k) to obtain
a(z-1)f(z-1)y(k)＝b(z-1)f(z-1)u(k)+w(k)；
Using the products obtained in step d
Figure BDA0001686893020000115
Calculating the input filtered signal:
Figure BDA0001686893020000116
and outputting a filtered signal:
Figure BDA0001686893020000117
and input and output filtered signals:
Figure BDA0001686893020000118
f. solving a secondary estimator of the parameters:
according to step d, finding
Figure BDA0001686893020000119
Then according to step e, obtain
Figure BDA00016868930200001110
And
Figure BDA00016868930200001111
by using
Figure BDA00016868930200001112
And
Figure BDA00016868930200001113
respectively in place of the formula f (z)-1)y(k)＝f(z-1)xT(k)θTY (k) and x in + w (k)T(k) And obtaining a common linear regression model between the output filtering signal and the input and output filtering signal:
Figure BDA00016868930200001114
when k is n +1, n +2, …, n + p, the vector matrix model is obtained:
Figure BDA00016868930200001115
wherein:
Figure BDA00016868930200001116
Figure BDA0001686893020000121
determining theta by using ordinary least square methodTSecond order approximation estimator of (1):
Figure BDA0001686893020000122
g. calculating the final parameter estimator:
by using
Figure BDA0001686893020000123
Substitution
Figure BDA0001686893020000124
According to step d, calculating
Figure BDA0001686893020000125
According to step f, calculating to obtain a third approximate estimator of the parameter
Figure BDA0001686893020000126
Repeating the iteration until
Figure BDA0001686893020000127
Wherein: i is
Figure BDA0001686893020000128
N is
Figure BDA0001686893020000129
The total number of iterations of (c);
at this time, θTIs estimated by
Figure BDA00016868930200001210
Then will be
Figure BDA00016868930200001211
Substitution into
Figure BDA00016868930200001212
And obtaining an estimated value or a predicted value of y (k), namely obtaining an estimated value or a predicted value of the forest accumulation amount of the k-th period of the target region.
From the above steps, the generalized least squares algorithm has the following characteristics: 1) in estimating ARMA model parameter vector thetaTMeanwhile, the noise model parameter vector f must be estimatedT(ii) a 2) Is a successive approximation algorithm, the specific calculation steps firstly assume fTAs is known, the ordinary least square method is used to find thetaTIs estimated by
Figure BDA00016868930200001213
Then use
Figure BDA00016868930200001214
Give fTIs estimated by
Figure BDA00016868930200001215
Reuse of
Figure BDA00016868930200001216
Correction removing
Figure BDA00016868930200001217
Repeating the iteration until
Figure BDA00016868930200001218
Until the end; 3) in each step of iteration process, the method is realized by adopting a common least square method.
The data of the forest land area and the forest accumulation amount in the peak town of Yongtai county, Fuzhou city in 1986-1999 are analyzed and explained, and the specific numerical values are shown in Table 1:
TABLE 1 urban town 1986-1999 forest area (1000 hm)2) Forest reserve (10000 m)3)
Figure BDA00016868930200001219
Figure BDA0001686893020000131
Results of calculations of one, different methods
Assuming that the time returns to 1997, the forest area and forest accumulation in 1986 and 1997 are the collected measured time series data, and the forest area in 1998 and 1999 is 0.523 and 0.534 ten thousand hm2Is predetermined by the township forest plan. The ARMA model was determined to be order 4. Thus, 4 methods were used to obtain estimates for the forest accumulation in the town of 1997 and 1990, whereas 4 methods were used to obtain predicted values for the forest accumulation in 1998 and 1999.
Let m be 2 (see equation (4) in step c), the common least squares parameter is:
(a1,a2,a3,a4,b0,b1,b2,b3,b4)＝(-0.3332,-0.1201,-0.5102,0.6819,1.2428,-0.0022,4.1678,0.2928,-3.1334)(f1,f2)＝(0.0006,-0.0024)
the generalized least squares parameter is:
(a1,a2,a3,a4,b0,b1,b2,b3,b4)＝(0.0816,-0.1260,-1.1647,0.1601,-2.8195,-1.3008,4.7038,2.8720,-3.7722)(f1,f2)＝(0,0)
the estimation and prediction results of the accumulation amount by using a feedforward BP neural network of a hidden layer are shown in a table 2, wherein the hidden layer comprises 5 neurons, and the output layer comprises 1 neuron.
Since the forest accumulation is a variable related to the area of the forest land, the gray model of GM (1,2) is used for estimation and prediction, and the GM (1,2) model is calculated as
Figure BDA0001686893020000132
Solving this differential equation yields the following function:
Figure BDA0001686893020000133
the actual value of the forest accumulation of the town, the estimated value and the predicted value of each method are shown in the table 2.
TABLE 2 evaluation of the accumulation and prediction results and comparison with the actual measurement values for the four methods
Figure BDA0001686893020000141
Second, result analysis
The residual variance of the noise parameter vector f, the estimated value and the predicted value obtained from the common least square algorithm can be obtained as follows: 1) the residual sequence xi (k) is colored noise, 2) the estimation of parameters by the ordinary least square is biased, so the estimation result of accumulation by the ARMA model and the ordinary least square is biased, and the prediction deviation of the forest accumulation in 1998 and 1999 is larger.
The residual variance of the noise parameter vector f, the estimated value and the predicted value obtained by the generalized least square algorithm can be obtained as follows: 1) the colored noise xi (k) has been converted into white noise w (k), 2) the estimation of parameters by the generalized least square is unbiased, so the estimation result of accumulation by the ARMA model and the generalized least square is unbiased, and the prediction precision of forest accumulation in 1998 and 1999 is the highest.
As can be seen from the estimation value of the GM (1,2) model, the residual variance of the predicted value and the prediction of the accumulation amount, the estimation and prediction effects of the GM (1,2) model are the worst, mainly because GM (1,2) is used
Figure BDA0001686893020000151
Fitting an exponential curve
Figure BDA0001686893020000152
Therefore, the estimation and prediction of the accumulation amount completely depend on the area of the wooded land. The accumulation is a time-related quantity, i.e. a certain autoregressive relation also exists between the accumulations in each periodThe grey model does not consider the autocorrelation among the accumulation quantities, and in addition, the common least square algorithm is adopted for parameter identification when the residual sequence is colored noise, so the estimation and prediction precision is the lowest.
The estimation value of the BP neural network, the residual variance of the predicted value, and the estimation and prediction of forest accumulation amount in 1990-: 1) BP neural network is very sensitive to initial weight and is very easy to converge on local minimum, 2) the determination of the number of the hidden nodes of the network has only some empirical formulas so far without any theoretical guidance.

Claims (1)

1. The unbiased prediction method of the forest accumulation amount is characterized by comprising the following steps: the method comprises the following steps:
1) collecting the paired data of the latest continuous p-period forest land area and forest accumulation amount actual measurement time sequence of the target region, wherein p is more than 2n +1, and n is the order of the model in the formula (1);
2) estimating the forest accumulation amount of the target region based on the autoregressive moving average model and the generalized least square algorithm through the data collected in the step 1), wherein the method specifically comprises the following steps:
a. a forest accumulation time sequence model and a noise lag operator polynomial:
the time series model of the forest accumulation is set as follows:
Figure FDA0001686893010000011
wherein: n is the order of the model;
k. j is constant and k is more than j, y (k) is the output of the system and represents the forest accumulation amount of the kth period;
y (k-j) is the historical output of the system and represents the forest accumulation amount in the k-j stage;
u (k-j) is a system determinability input and represents the forest land area at the k-j stage;
ajis an autoregressive coefficient, bjFor input of transfer coefficient, i.e. ajAnd bjParameters of a forest accumulation time series model;
ξ (k) is the residual error of the kth stage, i.e. the colored noise of the model;
let the noise lag operator polynomial be:
Figure FDA0001686893010000012
wherein: z is a radical of-i(i ═ 1, 2.. times.n) is a hysteresis operator, a1,a2,…,an,b0,b1,b2,…,bnIs a parameter;
equation (1) can be expressed as:
a(z-1)y(k)＝b(z-1)u(k)+ξ(k) (3)
b. determining the order n of the model:
let θ be (a)1,a2,...,an,b0,b1...,bn) Representing model parameter vectors; performing least square fitting by using models of different orders to calculate goodness of fit
Figure FDA0001686893010000021
Wherein
Figure FDA0001686893010000022
Is the parameter θ of the model at a certain orderTEstimating quantity;
when the order of the model is increased, the model is,
Figure FDA0001686893010000023
the value of (a) is reduced, and when the significant reduction is terminated, the order of the corresponding model is the appropriate model order;
c. calculating a first estimate of the parameter, i.e. the quantity that the generalized least squares algorithm needs to solve first:
let the whitening filter, i.e. the whitening function, be:
f(z-1)＝1+f1z-1+…+fmz-m； (4)
wherein: z is a radical of-j(j ═ 1, 2.. times, m) is a hysteresis operator, f1,f2,…,fmIs a parameter;
the conversion relation between the colored noise xi (k) and the white noise w (k) is
Figure FDA0001686893010000024
Wherein f (z)-1) Is a whitening filter;
let f (z) in formula (5)-1) 1, i.e. the noise model parameter vector f is an m-dimensional zero vector, where f ═ f1,f2,...,fm) Then, the residual error of the model formula (1) or formula (3) is degraded into white noise, i.e., ξ (k) ═ w (k);
order to
Figure FDA0001686893010000025
Wherein:
Figure FDA0001686893010000026
outputting information for the n-period history;
Figure FDA0001686893010000027
inputting information for the n +1 period;
when p pairs of sampling values are obtained
Figure FDA0001686893010000028
N +1, n +2, n + p, depending on y (k) xT(k)θT+ ξ (k), the vector matrix model is derived:
Y＝XθT+ξ； (6)
wherein:
YT＝(y(n+1),y(n+2),...,y(n+p))，
ξT＝(ξ(n+1),ξ(n+2),...,ξ(n+p))，
Figure FDA0001686893010000031
determining theta by ordinary least square methodTFirst approximation estimator of (1):
Figure FDA0001686893010000032
d. and (3) solving a noise model parameter vector:
order to
Figure FDA0001686893010000033
φT(k) Representing a colored noise sequence under random interference of a forest ecosystem; by using
Figure FDA0001686893010000034
And paired samples of each period, consisting of y (k) xT(k)θT+ ξ (k) the estimated value of colored noise ξ (k) at each stage of the model is calculated:
Figure FDA0001686893010000035
will obtain
Figure FDA0001686893010000036
Respectively substitute for
Figure FDA0001686893010000037
Obtaining:
Figure FDA0001686893010000038
according to the vector matrix model:
Figure FDA0001686893010000039
wherein:
Figure FDA00016868930100000310
Figure FDA00016868930100000311
wT＝(w(n+1),w(n+2),...,w(n+p))；
using the common least square method to obtain fTThe estimated amount of (c):
Figure FDA00016868930100000312
wherein:
Figure FDA00016868930100000313
is composed of
Figure FDA00016868930100000314
Transposing;
e. calculating a filtered signal:
substituting formula (5) into a (z)-1)y(k)＝b(z-1) u (k) + xi (k) to obtain
a(z-1)f(z-1)y(k)＝b(z-1)f(z-1)u(k)+w(k)；
Using the products obtained in step d
Figure FDA0001686893010000041
Calculating the input filtered signal:
Figure FDA0001686893010000042
and outputting a filtered signal:
Figure FDA0001686893010000043
and input and output filtered signals:
Figure FDA0001686893010000044
f. solving a secondary estimator of the parameters:
according to step d, finding
Figure FDA0001686893010000045
Then according to step e, obtain
Figure FDA0001686893010000046
And
Figure FDA0001686893010000047
by using
Figure FDA0001686893010000048
And
Figure FDA0001686893010000049
respectively in place of the formula f (z)-1)y(k)＝f(z-1)xT(k)θTY (k) and x in + w (k)T(k) And obtaining a common linear regression model between the output filtering signal and the input and output filtering signal:
Figure FDA00016868930100000410
when k is n +1, n +2, …, n + p, the vector matrix model is obtained:
Figure FDA00016868930100000411
wherein:
Figure FDA00016868930100000412
Figure FDA00016868930100000413
Figure FDA00016868930100000414
determining theta by using ordinary least square methodTSecond order approximation estimator of (1):
Figure FDA0001686893010000051
g. calculating the final parameter estimator:
by using
Figure FDA0001686893010000052
Substitution
Figure FDA0001686893010000053
According to step d, calculating
Figure FDA0001686893010000054
According to step f, calculating to obtain a third approximate estimator of the parameter
Figure FDA0001686893010000055
Repeating the iteration until
Figure FDA0001686893010000056
Wherein: i is
Figure FDA0001686893010000057
N is
Figure FDA0001686893010000058
The total number of iterations of (c);
at this time, θTIs estimated by
Figure FDA0001686893010000059
Then will be
Figure FDA00016868930100000510
Substitution into
Figure FDA00016868930100000511
And obtaining an estimated value or a predicted value of y (k), namely obtaining an estimated value or a predicted value of the forest accumulation amount of the k-th period of the target region.
CN201810575883.7A 2018-06-06 2018-06-06 Unbiased prediction method for forest accumulation Active CN108764583B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810575883.7A CN108764583B (en) 2018-06-06 2018-06-06 Unbiased prediction method for forest accumulation

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810575883.7A CN108764583B (en) 2018-06-06 2018-06-06 Unbiased prediction method for forest accumulation

Publications (2)

Publication Number Publication Date
CN108764583A CN108764583A (en) 2018-11-06
CN108764583B true CN108764583B (en) 2021-09-28

Family

ID=63999160

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810575883.7A Active CN108764583B (en) 2018-06-06 2018-06-06 Unbiased prediction method for forest accumulation

Country Status (1)

Country Link
CN (1) CN108764583B (en)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110795684A (en) * 2019-11-11 2020-02-14 南京喂啊游通信科技有限公司 Estimation method for observing noise variance
CN117494879B (en) * 2023-10-31 2024-06-14 四川大学 Aircraft taxiing time prediction method, device, equipment and readable storage medium

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103162676A (en) * 2013-04-10 2013-06-19 北京林业大学 Accurate estimation method for measuring forest sub compartment accumulation and biomass by utilizing angle gauge
CN103616014A (en) * 2013-11-29 2014-03-05 浙江农林大学 Method for measuring and computing forest stock volume parameter
CN105893737A (en) * 2016-03-24 2016-08-24 北京林业大学 Construction method of geological model-based forest growing stock estimation model
CN107909260A (en) * 2017-11-10 2018-04-13 浙江省地理信息中心 A kind of natural resource assets are left one's post audit commenting method

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8712897B2 (en) * 2012-05-31 2014-04-29 Hwey-Chyi LEE Stock analysis method, computer program product, and computer-readable recording medium
US10453026B2 (en) * 2015-03-06 2019-10-22 Walmart Apollo, Llc System and method for forecasting high-sellers using multivariate bayesian time series

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103162676A (en) * 2013-04-10 2013-06-19 北京林业大学 Accurate estimation method for measuring forest sub compartment accumulation and biomass by utilizing angle gauge
CN103616014A (en) * 2013-11-29 2014-03-05 浙江农林大学 Method for measuring and computing forest stock volume parameter
CN105893737A (en) * 2016-03-24 2016-08-24 北京林业大学 Construction method of geological model-based forest growing stock estimation model
CN107909260A (en) * 2017-11-10 2018-04-13 浙江省地理信息中心 A kind of natural resource assets are left one's post audit commenting method

Non-Patent Citations (10)

* Cited by examiner, † Cited by third party
Title
Allometric models to estimate biomass organic carbon stock in forest vegetation;Mohammed Alamgir 等;《Forestry Research》;20080706;第19卷(第2期);第101–106页 *
Generic allometric models including height best estimate forest biomass and carbon stocks in Indonesia;Ervan Rutishausera 等;《Forest Ecology and Management》;20131130;第307卷(第1期);第219-225页 *
Temporal Change in Aboveground Culms Carbon Stocks in the Moso Bamboo Forests and Its Driving Factors in Zhejiang Province;Xu,Lin 等;《Forests》;20171231;第8卷(第10期);第37页 *
一种面向下垫面不均一的森林碳通量监测方法;徐小军 等;《中国科学:信息科学》;20131015;第43卷(第10期);第1342-1352页 *
不同尺度下福建省杉木碳计量模型、预估及应用研究;林卓;《中国博士学位论文全文数据库 农业科技辑》;20160915;D049-3 *
基于BP神经网络模型的森林资源蓄积量动态估测;尹瑞安 等;《现代农业科技》;20180530;第154-155+161页 *
基于森林资源清查结果估算我国森林蓄积量变化;刘硕;《西北农林科技大学学报(自然科学版)》;20121031;第40卷(第11期);第147-151页 *
基于混合模型的油松林分蓄积量预测模型的建立;王少杰 等;《西北农林科技大学学报(自然科学版)》;20171208;第46卷(第02期);第147-151页 *
基于碳密度-林龄关系的黑龙江省森林碳汇潜力预测;任继勤 等;《环境科学研究》;20170415;第30卷(第04期);第552-558页 *
森林土壤碳库研究方法进展;周国模 等;《浙江林学院学报》;20060430;第23卷(第2期);第207-216页 *

Also Published As

Publication number Publication date
CN108764583A (en) 2018-11-06

Similar Documents

Publication Publication Date Title
Lee et al. Sparse markov decision processes with causal sparse tsallis entropy regularization for reinforcement learning
CN109002686B (en) Multi-grade chemical process soft measurement modeling method capable of automatically generating samples
Cherkassky et al. Computational intelligence in earth sciences and environmental applications: Issues and challenges
CN111967688B (en) Power load prediction method based on Kalman filter and convolutional neural network
CN110689179A (en) Water bloom prediction method based on space-time sequence mixed model
CN110309609B (en) Building indoor air quality evaluation method based on rough set and WNN
CN108764583B (en) Unbiased prediction method for forest accumulation
CN107122869A (en) The analysis method and device of Network Situation
CN111754034A (en) Time sequence prediction method based on chaos optimization neural network model
Hadi et al. The capacity of the hybridizing wavelet transformation approach with data-driven models for modeling monthly-scale streamflow
CN116628582A (en) Method for predicting concentration of gas in transformer oil based on online ensemble learning
CN115421216A (en) STL-ARIMA-NAR mixed model-based medium-and-long-term monthly rainfall forecasting method
CN114117852B (en) Regional heat load rolling prediction method based on finite difference working domain division
CN108776850B (en) Accurate estimation method for forest accumulation
Lola et al. Improving the performance of ann-arima models for predicting water quality in the offshore area of kuala terengganu, terengganu, malaysia
CN114970745B (en) Intelligent security and environment big data system of Internet of things
CN117131654A (en) Target observation method based on nonlinear optimal disturbance of pre-analysis initial guess condition
CN108827905B (en) near-infrared model online updating method based on local weighting L asso
CN115062764B (en) Intelligent illuminance adjustment and environmental parameter Internet of things big data system
CN112966449B (en) Manufacturing process quality defect root cause diagnosis method based on variable selection and probabilistic neural network
CN116307206A (en) Natural gas flow prediction method based on segmented graph convolution and time attention mechanism
Wang et al. Failure rate prediction based on AR model and residual correction
CN114580151A (en) Water demand prediction method based on gray linear regression-Markov chain model
CN114595890A (en) Ship spare part demand prediction method and system based on BP-SVR combined model
Ren et al. Short-term demand forecasting for distributed water supply networks: A multi-scale approach

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant