Estimation and bias correction of the magnitude of an abrupt level shift

(1)

Master Thesis in Statistics, Data Analysis and Knowledge

Discovery

Estimation and bias correction of the

magnitude of an abrupt level shift

(2)

(3)

Linköping University Electronic Press

Upphovsrätt

Detta dokument hålls tillgängligt på Internet – eller dess framtida ersättare –från

publiceringsdatum under förutsättning att inga extraordinära omständigheter

uppstår.

Tillgång till dokumentet innebär tillstånd för var och en att läsa, ladda ner,

skriva ut enstaka kopior för enskilt bruk och att använda det oförändrat för

icke-kommersiell forskning och för undervisning. Överföring av upphovsrätten vid en

senare tidpunkt kan inte upphäva detta tillstånd. All annan användning av

dokumentet kräver upphovsmannens medgivande. För att garantera äktheten,

säkerheten och tillgängligheten finns lösningar av teknisk och administrativ art.

Upphovsmannens ideella rätt innefattar rätt att bli nämnd som upphovsman i

den omfattning som god sed kräver vid användning av dokumentet på ovan

be-skrivna sätt samt skydd mot att dokumentet ändras eller presenteras i sådan form

eller i sådant sammanhang som är kränkande för upphovsmannens litterära eller

konstnärliga anseende eller egenart.

För ytterligare information om Linköping University Electronic Press se

för-lagets hemsida

http://www.ep.liu.se/

Copyright

The publishers will keep this document online on the Internet – or its possible

replacement –from the date of publication barring exceptional circumstances.

The online availability of the document implies permanent permission for

anyone to read, to download, or to print out single copies for his/hers own use and

to use it unchanged for non-commercial research and educational purpose.

Subsequent transfers of copyright cannot revoke this permission. All other uses of

the document are conditional upon the consent of the copyright owner. The

publisher has taken technical and administrative measures to assure authenticity,

security and accessibility.

According to intellectual property law the author has the right to be mentioned

when his/her work is accessed as described above and to be protected against

infringement.

For additional information about the Linköping University Electronic Press

and its procedures for publication and for assurance of document integrity, please

refer to its www home page:

http://www.ep.liu.se/.

(4)

Abstract

Consider a time series model which is stationary apart from a single shift in mean. If the time of a level shift is known, the least squares estimator of the magnitude of this level shift is a minimum variance unbiased estimator. If the time is unknown, however, this estimator is biased. Here, we first carry out extensive simulation studies to determine the relationship between the bias and three parameters of our time series model: the true magnitude of the level shift, the true time point and the autocorrelation of adjacent observations. Thereafter, we use two generalized additive models to generalize the simulation results. Finally, we examine to what extent the bias can be reduced by multiplying the least squares estimator with a shrinkage factor. Our results showed that the bias of the estimated magnitude of the level shift can be reduced when the level shift does not occur close to the beginning or end of the time series. However, it was not possible to simultaneously reduce the bias for all possible time points and magnitudes of the level shift.

Keywords: Change-point Detection, Bias Correction, Shrinkage Factor, Level Shift, Generalized Additive Model.

(5)

(6)

Acknowledgements

I would like to express my deep gratitude to my supervisor Prof. Anders Grimvall and Anders Norgaard who spent a lot of time on helping me develop the ideas and improve this thesis.

I also would like to thank my dear parents and all my good friends for their support and encouragement.

(7)

(8)

1 Introduction

1.1 Background

An abrupt change in the distribution of the observations in a time series is usually called a change-point. The most common type of change-point is a sudden level shift that is

superimposed on an otherwise constant or smoothly varying mean function. However, changes in dispersion of the collected data or in the slope of the mean function are also called change-points.

The origin of change point detection comes from the manufacturing industry where control plots were introduced as a basic detection method (Shewhart, 1924). Thereafter, change-point detection has spread to several other fields, such as bioinformatics, environmental science and the financial sector (Braun & Braun & Muller 2000; Alexandersson, 1986; Chen & Gupta, 1997). The methods that have been developed include numerical algorithms such as dynamic programming algorithms (Hawkins, 1977; Hawkins, 2001), likelihood ratio tests (Sen & Srivastava, 1975) and Bayesian methods for change-point detection (Chernoff & Zacks, 1964; Yao, 1984; Erdman & Emerson, 2008).

If the time of a level shift is known, the least squares estimator of the magnitude of this level shift is a minimum variance unbiased estimator. If the time is unknown, however, this estimator is biased. This has a simple explanation. If the true level shift occurs somewhere in the middle of the examined time interval, the least squares algorithm, in principle, searches for the time point that maximizes the difference in mean before and after, and this maximum is always greater than or equal to the difference in mean before and after the true change point. Changes in the mean that occur near the beginning or end of the time series are usually underestimated, because the least squares algorithm then tends to find a smaller but more centrally located shift in the observed data. It has previously been shown that the magnitude of a change detected by the CUSUM procedure can be strongly biased (Wu, 2004), but less attention has been paid to the bias of the ordinary less squares algorithm of a single level shift at an unknown time point.

Many statistical algorithms can be applied to correct biased estimates, such as bias-corrected jackknife and bootstrap estimators (Quenouille, 1956; Efron & Tibshirani, 1993) and the

(10)

2

penalization of the likelihood function (Firth, 1997). In this thesis, we first used the traditional plug-in method in which the bias is expressed as a function of some parameters and the true parameters are substituted for estimated parameters (MacKinnon & Smith, 1998). However, the finite sample properties can vary strongly with the shape of the bias function and the accuracy of the initial parameter estimates (MacKinnon & Smith, 1998). Therefore, we also introduced a shrinkage factor to reduce the bias.

1.2 Objectives

(1) To examine how the bias of the estimated magnitude of a level shift is influenced by the parameters of the investigated time series model.

(2) To investigate whether the bias can be reduced by modifying the standard least squares estimator.

(11)

3

2 Methodology

2.1 Least-squares estimation of a change point

The least-squares algorithm for estimating the time and magnitude of a single level shift splits a sequence of observations Y Y1, 2, ,Yn into two segments and computes the sum of squared

residuals (SSR) as: 2 2 1 2 1 1 ( ) ( ) ( ) t n i i i i t SSR t Y Y Y Y    



 



 where 1 1 Y t i i Y t  



and 1 2 Y n i i t Y n t    



The change point is defined as:

1

ˆ arg min (

( ))

t n

t

SSR t

 



and the estimated level shift is:

ˆ ( 1 ) ( 1 ˆ) ( 1 ˆ) ˆ ˆ ( ) n t t y y y y y y n t t            2.2 Bias Correction

The so-called plug-in method for bias correction (MacKinnon & Smith, 1998) based on the assumption that the expected value of a parameter estimate ˆ exists and can be expressed as:

ˆ

( )

( , ) /

E



 



b



γ

n

where n is the number of observations and b( , ) γ is a bias function in which γis a vector containing model parameters other than .

The corrected estimated parameter is

(12)

4

Besides the plug-in technique, a shrinkage factor was introduced to shrink the estimated parameter ˆ by computing

ˆ_Corr  ˆ ˆ* (2) where the shrinkage factor ˆwas calculated by minimizing the formula:

2 2 2

ˆ ˆ ˆ

( ) ( ) var

E_     _ E_{ }      _{ }

Formula 1 and 2 are the two basic methods for bias-correction used in this thesis.

2.3 Generalized additive models

The Generalized additive model was applied in this thesis to extrapolate the simulation results, the main reason for using GAM models is that they can provide the ability to detect the

nonlinear patterns without sacrificing interpretability (Hastie & Tibshirani, 1990). The Generalized Additive Model (GAM) was developed by Hastie and Tibshirani (1990).In this type of models, a function of the mean of the dependent variable is expressed as:

0 1 1 2 2

( ( )) ( ) ( ) _p( _p)

g E Y  s s x s x  s x

where g() is the link function, x₁, x_pare predictor variables and s₁( ), s_p( ) are unspecified smoothing functions. In practice s_p is estimated from the data by using techniques developed for scatterplot smoothing. There are various types of smoothers like running means, cubic splines, and B-splines. In this thesis, a combination of backfitting and B-spline smoothers were used to fit the model to data and the degrees of freedom of the smoothers was chosen by Generalized Cross Validation (Dong X. 2001; Wahba G. 1990).

(13)

5

3 Simulation experiments

3.1 Artificial datasets

To investigate how the bias is influenced by the true magnitude of the level shift, the location of the shift and the autocorrelation of the time series, artificial datasets were generated. More specifically, 75 groups of data were generated, each containing 10,000 time series with 100 values.

A single artificial time series:

A single time series data can be generated by following the steps below:

a.

Compute the constant _b ₁2

where  is the autocorrelation of adjacent observations.

b.

Lety₁ have a standard normal distribution and compute:

yi yi1   b oi



1, 100



1 (0,1) i (0,1)

i y N o N

To simulate series with one discontinuity, a step (level shift)  is artificially introduced in the series and its position is set as t₀, where t₀ has a uniform distribution on [1, 100].

(14)

6 100 80 60 40 20 0 5 4 3 2 1 0 -1 -2 -3 t y

Figure 1. Time series plot of 100 simulated values of a time series with a level shift of magnitude 3 at time 50. The autocorrelation at lag 1 was set to 0.3.

75 groups of artificial datasets:

A total of 75 groups (set the group number as g) of time series were computed by varying the magnitude and time of the change point and the autocorrelation coefficient (see Appendix). The magnitudes of the artificial level shift were set to 0.25, 0.5, 1, 2, and 3, and their positions (t₀)

were set to 10, 30, 50, 70 and 90. The artificial autocorrelation coefficients were set to 0, 0.3 and 0.6. Using this approach, 10,000 series of 100 values were generated for each combination of magnitude, position and autocorrelation, resulting in a total of 750,000 series with a single step.

3.2 Detection of change points in artificial datasets

For time series k within group g of the 75 groups described in the previous section, the

least-squares method was used to estimate the change-point

ˆt

as:

1

ˆ

_{arg min(}

_{( ))}

gk t n

t

SSR t

 



(15)

7 ˆ ˆ 1 1 1 ( ) ( ) ( ) ˆ ˆ ˆ ( ) gk gk n t t gk gk gk y y y y y y n t t           



1 75



g

Before estimating the autocorrelation coefficient _gk and the variance _gk of the series Y_gk the latter was homogenized by subtracting ˆ_gk after time point ˆt . Then, the variance was _gk

estimated as: 2 ( _i ) gk y y n   

and the first-order autocorrelation coefficient was set estimated as

1 (1) 1 (2) 1 1 1 1 ₂ 2 1 ₂ 2 (1) (2) 1 1 ( )( ) ˆ ( ) ( ) N i i i gk N N i i i i y y y y y y y y             _   _     



where y₍₁₎is the mean of the first n-1 observations from the homogenized series and y₍₂₎is the mean of the last n-1 observations from homogenized series, with n being the total number of points in the data series.

This process was run for each simulated time series, providing a total of 750,000 estimates. After sorting the results by group the average values ˆ_g, ˆt and g ˆg were computed

(16)

8

4 Bias Correction

4.1 Analysis of simulation results

The estimated bias of the magnitude of the level shift was calculated asBias  ˆ . Figures 2 and 3 illustrate the bias and dispersion of uncorrected estimates of the magnitude of level shifts. In figure 2, the true level shift occurred at time point 50 and had magnitude 0.5, and in figure 3 the magnitude of the level shift was 3 and it occurred at time point 10.

100 80 60 40 20 0 4 3 2 1 0 -1 -2 -3 -4

Estimated Change Point

Es ti m a te d L e v e l S h if t

Figure 2. 10 000 uncorrected estimates of the location and magnitude of a single level shift at time 50 in a time series of length 100. The true magnitude of the level shift was 0.5, as the solid line shows, and the mean value of the estimated level shift was 0.659 as the dashed line shows.

(17)

9 100 80 60 40 20 0 6 5 4 3 2 1 0 -1 -2 -3

Estimated Change Point

Es ti m a te d L e v e l S h if t

Figure 3. 10 000 uncorrected estimates of the location and magnitude of a single level shifts at time 10 in a time series of length 100. The true magnitude of the level shift was 3, as the solid line shows, and the mean value of estimated level shift was 2.910, as the dashed line shows.

From the plots in Figurs 2 and 3 we can see that the magnitude estimators are biased in both situations. These figures show that when the true level shift occurred in the middle of the time series, the magnitude of the level shift was slightly overestimated, whereas level shifts near the beginning or end of the time series were normally underestimated. In figure 2 the bias is

positive when the estimated time point of the level is in the interior of the time series, while the bias is more symmetric when the estimated time point is near the beginning or the end. In figure 3 the negative bias was obtained when the estimated time point is in the interior, while the estimates are more symmetrically spread at the end-points. So we can find that the bias does not only depend on where the true level shift is, but also on where the level shift is estimated to be and the bias of magnitude estimates could be a complex function of the true magnitude and location of the change point.

To illustrate the relationship between the bias and the true magnitude of the level shift, change point and data autocorrelation, the 750,000 bias estimates were grouped and plotted with different parameters:

(18)

10

Table 1. Descriptive statistics of the estimated bias of the magnitude of the level shift, grouped by the true level shift.

True level shift No. groups Mean Std Error Variance Minimum Median Maximum

0.25 15 0.0833 0.0148 0.0033 -0.0165 0.1069 0.1472

0.50 15 0.1183 0.0209 0.0065 -0.0462 0.1592 0.2051

1.00 15 0.0894 0.0173 0.0045 -0.0615 0.1021 0.1862

2.00 15 0.02286 0.00507 0.00039 -0.02631 0.02478 0.05095

3.00 15 0.00502 0.00176 0.00005 -0.00898 0.00552 0.01491

Table 2. Descriptive statistics of the estimated bias of the magnitude of the level shift, grouped by the change-point.

Change-point No. groups Mean Std Error Variance Minimum Median Maximum

10 15 0.0173 0.0114 0.0020 -0.0615 0.0149 0.1021

30 15 0.0910 0.0173 0.0045 -0.0062 0.1069 0.1826

50 15 0.0950 0.0199 0.0060 -0.0041 0.1079 0.2051

70 15 0.0851 0.0177 0.0047 -0.0016 0.1007 0.1792

90 15 0.0305 0.0104 0.0016 -0.0156 0.0257 0.1044

Table 3. Descriptive statistics of the estimated bias of the magnitude of the level shift, grouped by autocorrelation coefficients.

Autocorrelation

coefficients No. groups Mean Std Error Variance Minimum Median Maximum

0.0 25 0.0669 0.0111 0.0031 0.0016 0.0540 0.1599

0.3 25 0.0688 0.0127 0.0040 0.0024 0.0386 0.2021

(19)

11 0,185 0,160 0,135 0,110 0,085 0,060 0,035 0,010 t0 A rt ifi ci al le ve l s hi ft 90 80 70 60 50 40 30 20 10 3,0 2,5 2,0 1,5 1,0 0,5 0,185 -0,040 -0,015 0,010 0,035 0,060 0,085 0,110 0,135 0,160 bias 75 -0,1 50 0,0 0,1 3 25 0,2 2 1 0 0 Bias t0

Artificial level shift

Figure 4. Contour plot and 3D surface plot of bias versus the true level shift and change point.

0,110

0,085 0,010

ρ 3,0 2,5 2,0 1,5 1,0 0,5 0,6 0,5 0,4 0,3 0,2 0,1 0,0 0,185 -0,040 -0,015 0,010 0,035 0,060 0,085 0,110 0,135 0,160 bias

(20)

12 0 -0,1 0,0 0,1 0,6 0,2 1 _0,4 2 _{3 0,0} 0,2 Artificial level shift

Bias

ρ

Figure 5. Contour plot and 3D surface plot of bias of magnitude estimates versus true level shift and autocorrelation coefficient.

0,110 0,085 0,060 0,060 t0 ρ 90 80 70 60 50 40 30 20 10 0,6 0,5 0,4 0,3 0,2 0,1 0,0 0,185 -0,040 -0,015 0,010 0,035 0,060 0,085 0,110 0,135 0,160 bias -0,1 0,0 0,6 0,1 0,2 0 ₂₅ 0,4 0,2 50 75 0,0 Bias ρ t0

Figure 6. Contour plot and 3D surface plot of bias of magnitude estimates versus change point and autocorrelation coefficient.

(21)

13

The descriptive statistics tables (Tables 1,2 and 3) and the plots (Figures 4, 5 and 6) show that when the change point appears near the middle point of a time series and the true level is relatively small the bias reaches its maximum. Moreover, the presence of autocorrelation can increase the bias slightly.

The conclusion can be drawn that for all the three variables: artificial level shift, change point and autocorrelation coefficient, there is a noticeable difference between the cases in which bias is big and in which it is small, which suggests that all three variables might be useful in

building the bias function.

From Table 1, 2 and 3 we can find that the absolute bias of magnitude estimates is more frequently overestimated than underestimated so along with the bias estimates the shrinkage factor was computed to shrink the estimated level shift. To illustrate the relationship between the shrinkage factor and the true magnitude of the level shift, change point and data

autocorrelation, shrinkage estimates were grouped and plotted with different parameters using the same kind of graphs as for the bias analysis (Figures 7, 8 and 9).

0.96014 0.91289 0.86563 0.81837 0.77111 0.72386 0.67660 t0 A rt ifi ci al le ve l s hi ft 90 80 70 60 50 40 30 20 10 3.0 2.5 2.0 1.5 1.0 0.5 1.05466 1.10192 0.62934 0.67660 0.72386 0.77111 0.81837 0.86563 0.91289 0.96014 1.00740 sita

(22)

14 0.60 3 0.75 0.90 1.05 0 2 25 ₁ 50 75 0 θ

t0 Artificial level shift

Figure 7. Contour plot and 3D surface plot of shrinkage factors of magnitude estimates versus change point and true level shift.

0.96014 0.91289 0.81837

ρ 3.0 2.5 2.0 1.5 1.0 0.5 0.6 0.5 0.4 0.3 0.2 0.1 0.0 1.05466 1.10192 0.62934 0.67660 0.72386 0.77111 0.81837 0.86563 0.91289 0.96014 1.00740 sita 0.0 0.2 0.4 0.0 0.5 1.0 3 ₂ _0.6 1 ₀ ro θ

(23)

15

Figure 8. contour plot and 3D surface plot of shrinkage factors of magnitude estimates versus level shift and autocorrelation coefficient.

0.86563 t0 ρ 90 80 70 60 50 40 30 20 10 0.6 0.5 0.4 0.3 0.2 0.1 0.0 1.05466 1.10192 0.62934 0.67660 0.72386 0.77111 0.81837 0.86563 0.91289 0.96014 1.00740 sita 0.6 0.60 0.75 0.90 1.05 0.4 0 ₂₅ 0.2 50 75 0.0 θ ρ t0

Figure 9. Contour plot and 3D surface plot of shrinkage factors of magnitude estimates versus change point and autocorrelation coefficient.

From the plots in Figures 7, 8 and 9 we can find that all three variables might be useful to estimate the shrinkage factor.

4.2 Interpolation/ extrapolation of simulation results

The simulation results showed that the bias of magnitude estimates is a complex function of the data autocorrelation coefficient, the true magnitude and the location of the change point. We used generalized additive models to explore the relationships between the bias of the estimates and the model parameters.

(24)

16

A stepwise selection strategy was applied with different combinations of predictors to see if there is significant difference in the goodness of fit or deviance. The lower the deviance, the better the model fits the data.

Table 4. Expression of GAM models.

Type of Model: Mathematical Form:

Parametric E Y X( | x)01 1x 2x2

Nonparametric E Y X( | x)01xs x( )

Semiparametric E Y X( | x)01 1x 2x2s x( )1

Additive E Y X( | x)01 1x 2x2s x1( )1 s x2( 2)

Thin-plate spline E Y X( | x)0s x x( ,1 2)

Table 4 shows how to specify various models for a dependent variable y and independent variables x , x₁, and x₂. s_i( ) i1, 2 are nonparametric smooth functions,_i( ) i1, 2are parametric functions.

Bias model:

Table 5. Deviance for different combinations of predictors in the bias model

Model: Deviance: 0 1 1 0 2 ( ) ( ) ( ) ( ) E Bias    s t s  0.147 0 1 0 1 2 ( ) ( ) ( ) ( ) E Bias   t s  s  0.225 0 1 1 2 0 ( ) ( ) ( ) ( ) E Bias    s  s t 0.147 0 1 0 2 3 ( ) ( ) ( ) ( ) E Bias  s t s  s  0.146 0 1 0 2 ( ) ( ) ( , ) E Bias  s t s   0.129 0 1 2 0 ( ) ( ) ( , ) E Bias  s  s  t 0.095 0 1 2 0 ( ) ( ) ( , ) E Bias  s  s  t 0.113

Table 5 shows that the function E Bias( )0s1( ) s2( , ) t0 could get the smallest deviance,

(25)

17

Table 6: Analysis of deviance of the final bias model:

Source DF Sum of Squares Chi-Square Pr>Chi-Square Spline( ) 1.00 0.001 0.624 0.430

Spline(,t₀) 9.78 0.246 162.022 <0.0001

For each smoothing effect (spline) in the final model, Table 5gives a Chi-Square test comparing the deviance between the full model and the model without this variable. In this case, the results indicate that the effect of  and t₀is highly significant, whereas the effect of

 is insignificant at 5% level. The scatter plot in Figure 10 of predicted Bias versus Bias shows the model has a fairly good predictive ability.

0,20 0,15 0,10 0,05 0,00 -0,05 -0,10 0,20 0,15 0,10 0,05 0,00 -0,05 -0,10 Bias P re d ic te d B ia s

Figure 10. Scatterplot of Predicted Bias versus Bias.

From the simulation part we can conclude that (the value of) the shrinkage factor depends on the data autocorrelation coefficient, the true magnitude of the level shift and the location of the change point. We used generalized additive models to explore the relationships between the

(26)

18

shrinkage factor and the model parameters. A shrinkage model taking the shrinkage factor as the dependent variable, and _, t0 and as independent variables was estimated.

Shrinkage model:

Table 7: Deviance of the shrinkage model for different combinations of predictors.

Model: Deviance 0 1 1 2 0 ( ) ( ) ( ) ( ) E    _ s  s t 0.446 0 1 0 1 2 ( ) ( ) ( ) ( ) E    t s  s _ 0.565 0 1 1 0 2 ( ) ( ) ( ) ( ) E     s t s _ 0.350 0 1 0 2 3 ( ) ( ) ( ) ( ) E   s t s _ s  0.349 0 1 0 2 ( ) ( ) ( , ) E   s t s _  0.341 0 1 2 0 ( ) ( ) ( , ) E   s  s _ t 0.122 0 1 2 0 ( ) ( ) ( , ) E   s _ s  t 0.298

Table 7 shows that the function E( ) ₀s₁( ) s₂(_, )t₀ could get the smallest deviance, and this model was then chosen as the final bias function.

Table 8: Analysis of deviance of the final shrinkage factor model

Source DF Sum of Squares Chi-Square Pr>Chi-Square Spline( ) 1.00 0.002 0.748 0.387

Spline(_t₀) 16.73441 1.047 472.767 <0.0001

For each smoothing effect in the final shrinkage model, Table 8 gives a Chi-Square test

comparing the deviance between the full model and the model without this variable. The results indicate that the effect of _and is highly significant, whereas the effect of is

(27)

19 1.0 0.8 0.6 0.4 0.2 0.0 1.0 0.8 0.6 0.4 0.2 0.0 Shrinkage factor P re d ic te d s h ri n ka ge f ac to r

Figure 11. Scatterplot of predicted shrinkage factor ˆ versus .

4.3 Generalizing the bias function to arbitrary inputs

The SCORE statement in the GAM procedure was used to interpolate the simulation result. The SCORE procedure multiplies values from two data sets, one containing generalized additive models coefficients and the other containing raw data to be scored using the coefficients from the first data set. The scoring table includes 5 variables: artificial level shift, the location of change, correlation coefficient, expected bias and expected shrinkage factors. The artificial level shift was set to values ranging from 0.25 to 3 by increments of 0.05, the artificial change point t₀was set to values from 10 to 90 by increments of 1, and the artificial autocorrelation coefficient was set to values from 0 to 0.6 by increments of 0.1.

(28)

20

5 Correction results

The bias of the corrected level shift estimate is: Biascorr ˆcorr . Figure 12 shows that the

mean absolute bias decreased by more than 80 percent when bias correction was used.

75 70 65 60 55 50 45 40 35 30 25 20 15 10 5 1 0.20 0.15 0.10 0.05 0.00 -0.05 -0.10 Group Number B ia s

Mean of Corrected bias=0,0082

Mean of crude bias=0,0638

Crude bias Corrected bias

Figure 12. Original bias and the bias after correcting the magnitude estimates by using the bias function vs. group. The plot shows the bias grouped by artificial level shift, change point and data autocorrelation coefficient.

Figure 13 indicates that the bias function correction is adequate for some locations of the level shift while less suitable for others. For instance, when the change point occurred at time point 10, the absolute value of the corrected bias is 0.025 compared to the crude bias which is just 0.018.

(29)

21 3.00 2.00 1.00 0.50 0.25 0.12 0.09 0.06 0.03 0.00

M e a n V a lu e o f b ia s Crude bias Corrected bias 90 70 50 30 10 0.10 0.08 0.06 0.04 0.02 0.00 -0.02 -0.04 Change point M e a n V a lu e o f b ia s Crude bias Corrected bias 0,6 0,3 0,0 0.07 0.06 0.05 0.04 0.03 0.02 0.01 0.00 Aotucorrelation coefficient M e a n V a lu e o f b ia s Crude bias Corrected bias

Figure 13. Line plots of Mean values of crude bias and corrected bias with categorical variables as: artificial level shift, change point and autocorrelation coefficient.

75 70 65 60 55 50 45 40 35 30 25 20 15 10 5 1 0.2 0.1 0.0 -0.1 -0.2 Group Number B ia s

Mean of crude bias=0,0638

Mean of corrected bias=-0.0125

Crude Bias Corrected Bias

(30)

22

Figure 14. Comparing two sets of bias estimates, one set consisting of bias estimates after correcting the level shift estimates by using a shrinkage factor which (filled squares connected with a dashed line) , and the other set consisting of crude bias estimates (bulls connected with a solid line).

Both sets of bias estimates were grouped by artificial level shift, change point and data autocorrelation coefficient. 3.00 2.00 1.00 0.50 0.25 0.10 0.05 0.00 -0.05 M e a n V a lu e o f b ia s Crude Bias Corrected Bias

Artificial level shit

90 70 50 30 10 0.100 0.075 0.050 0.025 0.000 Change point M e a n V a lu e o f b ia s Crude bias Corrected bias 0.6 0.3 0.0 0.08 0.06 0.04 0.02 0.00 -0.02 -0.04 Autocorrelation coefficient M e a n V a lu e o f b ia s Crude bias Corrected Bias

Figure 15. Line plots of Mean values of crude bias and bias of shrink level shift estimator with categorical variables as: artificial level shift, change point and autocorrelation coefficient

From Figures 14 and 15 we can see that in general the shrinkage factor shrunk the level shift properly.

Table 9: Descriptive Statistics of original Bias, corrected bias A from bias function and corrected bias B from shrinkage model:

Number Mean SE Mean Variance Minimum Median Maximum Bias 75 0.06380 0.00787 0.00464 -0.06147 0.04260 0.20513

Bias A 75 0.00816 0.00611 0.00280 -0.11422 0.00070 0.12875

(31)

23

Comparing the mean value of the absolute bias we can see that the bias function normally produced better corrected results than the shrinkage factor

(32)

24

6 Discussion and Conclusion

This thesis started with a simulation study of the bias of the least squares estimator of the magnitude of a single level shift occurring at an unknown time point. The simulations showed that the bias is primarily influenced by the time and true magnitude of the level shift, and that the presence of autocorrelation can further increase the bias. They also showed that the

magnitude of the level shift can be strongly overestimated when the shift occurs not to close to the beginning or end of the time series.

In an attempt to remove the bias of magnitude estimates, a traditional plug-in method: bias function was applied. Another method to correct the bias is to introduce a shrinkage factor which is multiplied with the crude estimator of the magnitude of the level shift. The reason for using a shrinkage factor is that the absolute bias of magnitude estimates is more frequently overestimated than underestimated (Tables 1, 2, and 3). Generalized Additive Models have been used both to build the bias function and the model for the shrinkage parameter. From the simulation results we have found that the variables (bias and shrinkage factor) have nonlinear relationships with the independent variables so it would not have been optimal to use

Generalized Linear Models. Other statistical techniques could have been applied, such as Neural Network, which have shown to work successfully for nonlinear modeling. However, the trade-off would be the lack of interpretability.

After have chosen the proper Generalized Additive Models, the score table was build to

interpolate the simulation results and the bias was corrected. From the correction results we can see that the bias function did not perform well in some situations since the variance of the estimates of the change points is very large (Figure 2). In general, there is a trade-off between bias and variance of parameter estimates (MacKinnon & Smith, 1998). However, our bias function and shrinkage factor both reduced the bias and the mean square error.

At last, three conclusions can be drawn. First, the bias of estimates is strongly influenced by the location of the level shift. Second, the estimation of strongly significant level shift is almost unbiased while less significant level shifts could lead to uncertain parameter estimates and biased results. Third, the simulations showed that the bias and mean square error can be substantially reduced by multiplying the level shift estimate by a shrinkage factor.

(33)

25

7 References

Alexandersson H. 1986. A homogeneity test applied to precipitation data. Journal of Climatology 6:661–675.

Braun J V, Braun R K and Muller H G. 2000. Multiple change-point ﬁtting via quasilikelihood, with application to DNA sequence segmentation. Biometrika 87:301–314.

Chen J and Gupta A K. 1997. Testing and locating variance change points with application to stock prices, Journal of American Statistical Association 92:739-747.

Dong X. 2001. Fitting Generalized Additive Models with the GAM procedure. SUGI Proceedings.

Erdman C and Emerson J W. 2008. A fast Bayesian change point analysis for the segmentation of microarray data. Bioinformatics 24:2143-2148.

Elmasri R and Navathe B. 2007. Fundamentals of Database Systems. Addison Wesley.

Firth D. 1997. Bias reduction of maximum likelihood estimates. Biometrika 80:27-38.

Hastie T and Tibshirani R. 1986. Generalized additive models (with discussion). Statistical Science 1:297-318.

Hastie T and Tibshirani R. 1987. Generalized additive models: some applications. Journal of American

Statistical Association 82:371-386.

Hastie T, Tisbshirani R and Friedman J. 2001. The elements of statistical learning: Data Mining,

Inference, and Prediction. New York: Springer Verlag.

Hawkins D M. 1977. Testing a sequence of observations for a shift in location. Journal of American

Statistical Association 72:180–186.

Hawkins D M. 2001. Fitting multiple change-point models to data. Computational Statistics and Data

Analysis 37:323–341.

Johnson A and Wichern W. 2003. Applied Multivariate Statistical Analysis. Chicago: Pearson Education.

Jushan B 1997. Estimation of a change point in multiple regression models. The Review of Economics

and Statistics 79:551-563.

Khaliq M N and Ouarda T B M J. 2007. On the critical values of the standard normal homogeneity test (SNHT). International Journal of Climatology 27:681-687.

Lora Delwiche and Susan Slaughter. 2008. The Little SAS Books: A Primer. SAS Publishing.

MacKinnon J G and Smith Jr A A. 1998. Approximate bias correction in econometrics. Econometrics 85:205-230.

Quenouille M H. 1956. Notes on bias in estimation. Biometrika 43:353-360.

(34)

26

Wahba G. 1990. Spline Models for Observational Data. Philadelphia: Society for Industrial and Applied Mathematics.

Wu Y. 2004. Bias of estimator of change point detected by a CUSUM procedure. Annals of The Institute

of Statistical Mathematics 56:127-142.

Yao Y. 1984. Estimation of a noisy discrete-time step function: Bayes and empirical Bayes approaches. The Annals of Statistics 12: 1434-1447.

(35)

27

8 Appendix

The detection results of 75 group simulation data.

Group

No. t0    ˆ bias

ˆ ( )

MSE

Correction results from bias function

Correction results from shrinkage factor

ˆ

corr

 ₍ˆ ₎

g

MSE ˆ_corr MSE(ˆ_g)

1 10 3.00 0.0 1.00 3.01 0.01 0.13 3.01 0.13 2.99 0.15 2 10 2.00 0.0 0.98 2.04 0.04 0.14 2.02 0.14 2.05 0.14 3 10 1.00 0.0 0.91 1.10 0.10 0.39 1.06 0.40 0.76 0.48 4 10 0.50 0.0 0.87 0.57 0.07 0.81 0.51 0.83 0.33 0.63 5 10 0.25 0.0 0.85 0.29 0.04 0.98 0.22 1.00 0.16 0.67 6 30 3.00 0.0 1.00 3.01 0.01 0.05 2.99 0.05 3.08 0.09 7 30 2.00 0.0 0.99 2.02 0.02 0.05 1.96 0.05 1.83 0.14 8 30 1.00 0.0 0.94 1.07 0.07 0.06 0.97 0.07 0.57 0.30 9 30 0.50 0.0 0.76 0.66 0.16 0.41 0.56 0.39 0.23 0.39 10 30 0.25 0.0 0.66 0.38 0.13 0.83 0.30 0.82 0.13 0.59 11 50 3.00 0.0 1.00 3.00 0.00 0.04 3.00 0.04 3.08 0.07 12 50 2.00 0.0 0.99 2.01 0.01 0.04 1.95 0.05 1.83 0.14 13 50 1.00 0.0 0.95 1.05 0.05 0.04 0.94 0.05 0.55 0.31 14 50 0.50 0.0 0.76 0.66 0.16 0.32 0.55 0.30 0.20 0.34 15 50 0.25 0.0 0.63 0.39 0.14 0.78 0.31 0.76 0.12 0.55 16 70 3.00 0.0 0.99 3.00 0.00 0.05 2.99 0.05 3.05 0.08 17 70 2.00 0.0 0.99 2.01 0.01 0.05 1.96 0.05 1.80 0.16 18 70 1.00 0.0 0.95 1.06 0.06 0.07 0.96 0.07 0.54 0.32 19 70 0.50 0.0 0.76 0.66 0.16 0.44 0.57 0.42 0.23 0.42 20 70 0.25 0.0 0.66 0.38 0.13 0.87 0.30 0.86 0.14 0.63 21 90 3.00 0.0 1.00 3.00 0.00 0.10 3.00 0.11 3.01 0.18 22 90 2.00 0.0 0.99 2.03 0.03 0.11 2.00 0.11 1.77 0.29 23 90 1.00 0.0 0.91 1.10 0.10 0.32 1.05 0.33 0.69 0.47 24 90 0.50 0.0 0.83 0.60 0.10 0.80 0.54 0.80 0.33 0.64 25 90 0.25 0.0 0.82 0.30 0.05 1.02 0.24 1.03 0.16 0.72 26 10 3.00 0.3 1.00 3.01 0.01 0.22 3.01 0.22 2.97 0.22 27 10 2.00 0.3 0.98 2.04 0.04 0.26 2.02 0.27 1.81 0.46 28 10 1.00 0.3 0.95 1.06 0.06 0.61 1.01 0.64 0.73 0.62 29 10 0.50 0.3 0.95 0.52 0.02 0.88 0.45 0.88 0.31 0.58 30 10 0.25 0.3 0.97 0.26 0.01 0.97 0.19 0.99 0.15 0.56 31 30 3.00 0.3 1.00 3.01 0.01 0.09 2.99 0.09 3.01 0.14 32 30 2.00 0.3 0.99 2.02 0.02 0.09 1.97 0.10 1.76 0.25 33 30 1.00 0.3 0.89 1.13 0.13 0.16 1.03 0.15 0.58 0.37 34 30 0.50 0.3 0.74 0.68 0.18 0.55 0.59 0.53 0.26 0.40

(36)

28 35 30 0.25 0.3 0.67 0.37 0.12 0.86 0.29 0.84 0.14 0.50 36 50 3.00 0.3 1.00 3.01 0.01 0.07 2.99 0.08 3.01 0.12 37 50 2.00 0.3 0.99 2.02 0.02 0.07 1.96 0.08 1.76 0.25 38 50 1.00 0.3 0.90 1.11 0.11 0.10 1.00 0.10 0.55 0.36 39 50 0.50 0.3 0.71 0.70 0.20 0.47 0.60 0.44 0.25 0.37 40 50 0.25 0.3 0.63 0.40 0.15 0.82 0.31 0.80 0.14 0.48 41 70 3.00 0.3 1.00 3.00 0.00 0.08 2.99 0.09 2.99 0.14 42 70 2.00 0.3 0.99 2.02 0.02 0.08 1.97 0.09 1.73 0.28 43 70 1.00 0.3 0.90 1.11 0.11 0.15 1.02 0.14 0.56 0.37 44 70 0.50 0.3 0.74 0.68 0.18 0.56 0.59 0.54 0.27 0.42 45 70 0.25 0.3 0.67 0.37 0.12 0.89 0.29 0.88 0.15 0.53 46 90 3.00 0.3 1.00 3.00 0.00 0.18 3.00 0.19 2.94 0.30 47 90 2.00 0.3 0.99 2.03 0.03 0.21 2.00 0.22 1.73 0.48 48 90 1.00 0.3 0.93 1.08 0.08 0.52 1.03 0.54 0.69 0.59 49 90 0.50 0.3 0.90 0.56 0.06 0.87 0.49 0.88 0.32 0.59 50 90 0.25 0.3 0.89 0.28 0.03 0.99 0.21 1.01 0.16 0.59 51 10 3.00 0.6 1.00 2.99 -0.01 0.43 2.99 0.44 2.76 0.58 52 10 2.00 0.6 1.01 1.97 -0.03 0.61 1.95 0.64 1.65 0.85 53 10 1.00 0.6 1.07 0.94 -0.06 0.99 0.89 1.05 0.64 0.78 54 10 0.50 0.6 1.10 0.45 -0.05 1.10 0.39 1.14 0.28 0.59 55 10 0.25 0.6 1.07 0.23 -0.02 1.13 0.16 1.15 0.14 0.53 56 30 3.00 0.6 0.98 3.01 0.01 0.18 3.00 0.19 2.89 0.28 57 30 2.00 0.6 0.98 2.05 0.05 0.19 2.00 0.18 1.69 0.45 58 30 1.00 0.6 0.85 1.18 0.18 0.41 1.11 0.40 0.65 0.45 59 30 0.50 0.6 0.75 0.67 0.17 0.85 0.59 0.83 0.31 0.46 60 30 0.25 0.6 0.70 0.36 0.11 1.05 0.28 1.05 0.16 0.49 61 50 3.00 0.6 1.00 3.01 0.01 0.15 3.00 0.15 2.88 0.25 62 50 2.00 0.6 0.98 2.04 0.04 0.14 2.00 0.15 1.68 0.41 63 50 1.00 0.6 0.84 1.19 0.19 0.30 1.11 0.29 0.63 0.41 64 50 0.50 0.6 0.71 0.71 0.21 0.08 0.63 0.75 0.31 0.43 65 50 0.25 0.6 0.63 0.38 0.13 1.05 0.31 1.03 0.17 0.50 66 70 3.00 0.6 1.00 3.01 0.01 0.17 3.00 0.18 2.86 0.28 67 70 2.00 0.6 0.98 2.04 0.04 0.17 2.00 0.19 1.67 0.45 68 70 1.00 0.6 0.85 1.17 0.17 0.38 1.10 0.37 0.64 0.45 69 70 0.50 0.6 0.75 0.67 0.17 0.84 0.60 0.83 0.31 0.47 70 70 0.25 0.6 0.71 0.35 0.10 1.06 0.28 1.05 0.16 0.50 71 90 3.00 0.6 1.00 2.99 -0.01 0.36 2.99 0.37 2.81 0.57 72 90 2.00 0.6 1.00 2.00 0.00 0.46 1.98 0.49 1.64 0.78 73 90 1.00 0.6 1.00 1.00 0.00 0.89 0.95 0.92 0.65 0.72 74 90 0.50 0.6 1.03 0.48 -0.02 1.09 0.42 1.12 0.30 0.59 75 90 0.25 0.6 1.01 0.25 0.00 1.13 0.18 1.15 0.15 0.53

(37)

Rapporttyp Report category Licentiatavhandling Examensarbete C-uppsats D-uppsats Övrig rapport _ ________________ Språk Language Svenska/Swedish Engelska/English _ ________________ Titel Title

Estimation and bias correction of the magnitude of an abrupt level shift

Författare Author

Wenjie Liu

Sammanfattning Abstract

Consider a time series model which is stationary apart from a single shift in mean. If the time of a level shift is known, the least squares estimator of the magnitude of this level shift is a minimum variance unbiased estimator. If the time is unknown, however, this estimator is biased. Here, we first carry out extensive simulation studies to determine the relationship between the bias and three parameters of our time series model: the true magnitude of the level shift, the true time point and the autocorrelation of adjacent observations. Thereafter, we use two generalized additive models to generalize the simulation results. Finally, we examine to what extent the bias can be reduced by multiplying the least squares estimator with a shrinkage factor. Our results showed that the bias of the estimated magnitude of the level shift can be reduced when the level shift does not occur close to the beginning or end of the time series. However, it was not possible to simultaneously reduce the bias for all possible time points and magnitudes of the level shift.

ISBN

_____________________________________________________ ISRN

LIU-IDA/STAT-A--12/007—SE

Serietitel och serienummer ISSN

Title of series, numbering ____________________________________

Nyckelord Keyword

Change-point Detection, Bias Correction, Shrinkage Factor, Level Shift, Generalized Additive Model. Datum

Date 10/14/2012

URL för elektronisk version

http://urn.kb.se/resolve?urn=urn:nbn:se: liu:diva-84618

Avdelning, Institution Division, Department

Division of Statistics, Department of Computer and Information Science

Division, Department

Division of Statistics, Department of Computer and Information Science

(38)

Estimation and bias correction of the magnitude of an abrupt level shift

Master Thesis in Statistics, Data Analysis and Knowledge

Discovery