• No results found

Wiener and Gamma Processes Overview for Degradation Modelling and Prognostic

N/A
N/A
Protected

Academic year: 2022

Share "Wiener and Gamma Processes Overview for Degradation Modelling and Prognostic"

Copied!
83
0
0

Laster.... (Se fulltekst nå)

Fulltekst

(1)

Wiener and Gamma Processes Overview for Degradation Modelling and Prognostic

Yun Zhang

Reliability, Availability, Maintainability and Safety (RAMS) Supervisor: Anne Barros, IPK

Department of Production and Quality Engineering Submission date: June 2015

Norwegian University of Science and Technology

(2)
(3)

RAMS

Reliability, Availability, Maintainability, and Safety

Wiener and Gamma Processes Overview for Degradation Modelling and Prognostic

Yun Zhang

June 2015

TPK4550 MASTER THESIS

Department of Production and Quality Engineering Norwegian University of Science and Technology

Supervisor: Professor Anne Barros

(4)

i

Preface

This master thesis was written during the spring semester 2015 at the Department of Production and Quality Engineering, Faculty of Engineering, Science and Technology, Norwegian University of Science and Technology. The thesis is part of the two-year international master’s program Reliability, Availability, Maintainability and Safety (RAMS).

The topic, which brought up in October 2014, was initally condition-based maintenance for stochastically deteiorating system. It was later slightly narrowed to become a study on Gamma process in degradation modeling and related maintenance policies. The title of this thesis is

"Wiener and Gamma Processes Overview for Degradation Modelling and Prognostic" and is written with the guidance of my supervisor professor Anne Barros at the Department of Pro- duction and Quality Engineering.

This report is partly based on the theory and methods in TPK4120 Safety and Reliability Anal- ysis, TMA4275 Lifetime Analysis and in TPK5170 RAMS Assessment and Optimization. So it is written for readers with basic knowledge in probability, reliability analysis, stochastic process theory, and preferably maintenance optimization.

Trondheim, 2015-06-08

Yun Zhang

(5)

ii

Acknowledgment

I would like to express my deepest gratitude to my supervisor, Professor Anne Barros, Depart- ment of Production and Quality Engineering at NTNU, for all her excellent guidance and sup- port throughout my master thesis. She provided me valuable advices and much of the informa- tion needed to write this thesis.

Y.Z.

(6)

iii

Summary

High reliability is an indispensable requirement for the operation of technical systems and in- frastructure (like roads, railways, buildings, bridges and industrial plants). Failures in these ar- eas can result in high costs and great hazards to humans and the environment. Many failure mechanisms can be traced to an underlying degradation process. Therefore, inspections and condition-based maintenance (CBM) are undertaken to monitor deterioration and to prevent damage and future failures.

Prognostic of component lifetime is important for CBM in many application domains where safety, reliability and availability are considered of first importance. However, one of the condi- tion indexes - the remaining useful lifetime (RUL) prediction is seldom taken into account in the decision making and maintenance planning in practice. To make up for it, the thesis focuses on a stochastic degradation process for RUL estimation and prognostic use.

The thesis starts with reviewing some of the degradation models. The merits, limitations of each model are presented. By aggregating the information of each model, this paper provides the key information about circumstances for choosing suitable deterioration models in the con- text of maintenance optimization.

Two stochastic process - Brownian motion and Gamma process are discussed in detail. Their statistical properties, methods for estimation, and simulation of are systematically reviewed with numerical examples. Also, each of the model is associated with the component’s uncer- tainties in the observations while estimating its RUL.

Since Gamma process has proven to be very useful in modeling degradation paths, deter- mining optimal inspection and maintenance decisions, the paper then investigates the appli- cation of Gamma degradation process in maintenance policies. Two RUL related maintenance policies are proposed and compared with a traditional degradation level-based policy. The per- formances of the proposed policies are evaluated through numerical examples in previous pa- pers and their advantages are stated in the end.

(7)

Contents

Preface . . . i

Acknowledgment . . . ii

Summary. . . iii

1 Introduction 1 1.1 Background . . . 1

1.2 Objectives . . . 4

1.3 Scope and Limitations . . . 5

1.4 Approach . . . 5

1.5 Structure of the Report . . . 5

2 Review of the Problem 7 2.1 Literature Review on Degradation Modeling and Prognostic Approaches . . . 7

2.2 Degradation Modeling . . . 10

2.2.1 Random coefficient regression Models . . . 10

2.2.2 Brownian Motion . . . 11

2.2.3 Gamma Process . . . 12

3 Brownian Motion Model 14 3.1 Standard Brownian Motion - Wiener Process . . . 14

3.1.1 Introduction and Definition. . . 14

3.1.2 Properties of Standard Brownian Motion . . . 17

3.1.3 Properties of Standard Brownian Path . . . 18

3.2 Brownian Motion with Linear Drift. . . 19

iv

(8)

CONTENTS v

3.2.1 Properties of Brownian Motion with Linear Drift . . . 19

3.2.2 Some Graphs and Discussions . . . 21

3.2.3 First Hitting Time and RUL Distribution . . . 23

3.3 Philosophy of RUL Estimation . . . 26

3.3.1 Assumptions and Training Data Set . . . 27

3.3.2 Parameter Estimation . . . 27

3.3.3 Likelihood Confidence Interval and Likelihood Test . . . 30

3.3.4 RUL Distribution . . . 32

3.4 Discussion about Model Relevance . . . 33

4 Gamma Process Model 35 4.1 Homogeneous Gamma Process . . . 35

4.1.1 Introduction and Definition. . . 35

4.1.2 Properties of Homogeneous Gamma Process . . . 36

4.1.3 Some Graphs and Discussions . . . 38

4.1.4 First Hitting time and RUL distribution . . . 41

4.2 Philosophy of RUL Estimation . . . 43

4.2.1 Assumptions and Training Data Set . . . 43

4.2.2 Parameter Estimation . . . 44

4.2.3 Likelihood Confidence Interval . . . 47

4.2.4 RUL Distribution . . . 47

4.3 Discussion about Model Relevance . . . 48

5 Condition-based Maintenance with Gamma-Process Deterioration 50 5.1 Assumptions and Performance Assessment Criterion. . . 51

5.2 Degradation-based Maintenance Model (T,M) Policy . . . 52

5.3 CDF RUL-based Maintenance Model (T,Q) Policy . . . 53

5.4 Mean RUL-based Maintenance Model (T,τm) Policy . . . 54

5.5 Comparison and Performance Assessment . . . 55

6 Summary 56 6.1 Summary and Conclusions . . . 56

(9)

CONTENTS vi

6.2 Discussion and Recommendations for Further Work . . . 57

A Acronyms 60 B Definitions and Theorems 61 B.1 Stopping Time . . . 61

B.2 Strong Markov Property . . . 61

B.3 New Process . . . 61

C Matlab Code for Wiener Process Simulation in Chapter 3 62 C.1 Example 1 . . . 62

C.2 Example 2 . . . 63

C.3 Example 3 . . . 63

C.4 Testing Data Generation and Parameter Estimation . . . 64

D Matlab Code for Gamma Process Simulation in Chapter 4 66 D.1 Realizations of gamma process path . . . 66

D.2 Testing Data Generation and Parameter Estimation . . . 67

Bibliography 70

(10)

Chapter 1 Introduction

This chapter presents a description of the background for the project thesis first. The thesis objectives are then defined, along with its scope, limitations and approach. Structure of the report is provided at the end of this chapter.

1.1 Background

Many technical systems and engineering infrastructures are subject to failures resulting from deterioration with usage and age. Power plants, railway tracks, oil platforms, cutting tools, bridges, dike and pipelines are some examples to illustrate such degradation processes (Jeang (1999)). The degradation and failures of these systems could incur high costs (e.g. due to down- time production losses and unplanned maintenance actions on the systems) and pose hazards to people and the environment. In order to reduce the costs, practices like inspections and maintenance models are increasingly applied in different fields.

An important feature of maintenance optimization is that decision often must be made un- der uncertainty, such as in degradation and maintenance cost (Van Noortwijk (2009)). If the evolution of a certain component or system degradation is explicitly known or even adequately predicted, optimum scheme and better decisions could be carried out to maintain the degraded system before failures. For this reason, a growing interest has been drawn to model the degra- dation under uncertainty for precise prognostics and wise maintenance measurements.

According toJardine et al.(2006), there are two main prediction types in prognostics. One is 1

(11)

CHAPTER 1. INTRODUCTION 2 to predict how much time is left before a failure occurs given the current subject condition and history operation profile. The length from the current time to the end of the useful life is called remaining useful life (RUL). Another type is to predict the chance that a subject operates without a fault or a failure up to some future time given the history and current status. Although there are more and more papers (see for exampleJardine et al.(2006),Welte et al.(2006),Van Noortwijk (2009), Liu et al.(2013))covering many aspects of prognostic problems on RUL over the last few decades, most methods exhibit limits in cases with uncertainty from the inner states or the external operating conditions of systems.

Zhang et al.(2015) points out that such kind of uncertainty (in his words heterogeneity) is wide spread in both the inner states of the system and the related working environments. For instance, an aircraft engine may experience various operational modes, such as lift, drag and take off etc. Also, the engine may perform differently under different workloads. Even com- ponents from the same category may exhibit various degradation paths in the same environ- ment. Therefore, in order to achieve a more accurate RUL estimation, we need to incorporate the above heterogeneity into degradation modeling. Towards this point, it is desirable to clas- sify heterogeneity into three categories: the unit-to-unit variability for components from the same category, the variability in time-varying operating conditions, and the diversity of tasks and workloads of systems during their life cycles. The taxonomy of RUL estimation approaches for system with heterogeneity is illustrated in Fig.1.1.

From the figure we can see that RUL and the associated degradation modeling techniques are composed of large branches of situations where their uncertainties are quite different from each other. Adding one kind of heterogeneity to another could increase the dimension of com- plexity in degradation modeling. However, the combination of heterogeneity in RUL estima- tion is relatively novel and unexplored. There are few comprehensive papers on combined ap- proaches for prognostics with RUL estimation. Therefore, the writer aims to improve the knowl- edge in each kind of heterogeneity and propose a hybrid solution for degradation modeling based RUL estimation for systems with three heterogeneity. As a starting point, this master the- sis focuses on the first kind of heterogeneity (i.e. methods considering unit-to-unit variability) and only studies RUL estimation from the statistical point of view (i.e. stochastic process with random effects).

(12)

CHAPTER 1. INTRODUCTION 3

Figure 1.1: Taxonomy of RUL estimation approaches for system under heterogeneity (adopted fromZhang et al.(2015) ).

The concept of RUL has been in various fields ranging from material science, biostatistics and econometrics (Si et al.(2011)). The word "useful" is normally economic related and de- fined upon individual explanations within specific context and operational characteristics. In this thesis, the writer assume that the definition of the useful life is know to the owner of the asset and the interest is to investigate the methods for degradation modeling related to RUL estimation given condition information.

Prognostic based on RUL estimation is one of the key factors in condition based main- tenance, and prognostics and health management (Cui et al. (2004), Lee et al. (2006), Wang (2007b), Wang (2007a)). The estimation results could have great impacts on the planning of maintenance activities, spare parts provision, operational performance and the profitability of the owner of an asset(Jardine et al.(2006),Altay and Green(2006),Elwany and Gebraeel(2008), Papakostas et al.(2010)). For example, one may need to estimate the remaining technical life time of an industrial machine without data on failures. In other cases, one may want to sched- ule future repairs and maintenance on when to replace the machine, etc. RUL estimation has also an important role in the reusing and recycling products which , as a result, has a strategic impact on energy consumption, raw material use, pollution, and landfill (Mazhar et al.(2007)).

(13)

CHAPTER 1. INTRODUCTION 4 Therefore, it is very important to estimate RUL beyond CBM and prognostics.

1.2 Objectives

The thesis focuses on statistical-based approaches for degradation modeling and prognostic considering only unit-to-unit variability. The approaches rely on available past monitoring data and stochastic models to estimate the RUL in a probabilistic way. Such method is based on the fact that condition monitoring data and extracted features vary with the development of either the initiation and propagation process or the degradation process (Ahmadzadeh and Lundberg (2014)). So it requires no special product knowledge, physics or engineering principles. The method also have certain advantages as some nice mathematical properties can be analyzed regarding to the estimated RUL. The investigation result is intended to serve as a basis for main- tenance optimization combining the merits of stochastic models (like Gamma process). This master thesis can be led with some connection with SUBPRO project. To achieve all these goals, the following tasks are identified:

1. Make a brief survey about degradation modeling to complement the work achieved dur- ing the specialized project.

2. Improve the knowledge to be able to implement the whole modeling process on some case studies published in previous works.

3. Specially focus on parameter estimation step when using degradation models consider- ing stochastic process with random effects (e.g. Gamma process), in case of perfect and imperfect observations (e.g. white noise).

4. Study the Remaining Useful Lifetime (i.e. calculation of its law and estimation of the rele- vant parameters) in case of monotone and non-monotone degradation models.

5. Transport the merits of degradation modeling to maintenance optimization.

(14)

CHAPTER 1. INTRODUCTION 5

1.3 Scope and Limitations

This thesis is directed towards students, professors and researchers and other people who have basic knowledge of RAMS engineering and carry out required analysis in the condition-based maintenance and prognostic. The thesis is mainly an investigation and also a summary on degradation modeling with a focus on stochastic process. It serves as a theoretical preparation for related maintenance optimization policy. Most parts of this thesis involves mathematical models and statistical properties, and the development of an optimal maintenance policy is not quantitatively discussed here due to limited amount of time. However, it should be empha- sized that the interest of this thesis is not only in mathematics, but more importantly, in applied statistics and RAMS optimization.

1.4 Approach

A large proportion of the thesis is established based on literature survey and text books. Simu- lation are implemented through out the thesis to demonstrate and illustrate the modeling pro- cesses, and it enhances the understanding of stochastic process and degradation modeling. Ob- jective 1 will be established through the literature review from many RUL related papers scat- tered among operational research, reliability modeling, optimal maintenance, fault diagnosis and prognosis, and survival analysis literature. Objective 2,3,4 will be approached according to books and articles about stochastic process, lifetime analysis and applied statistics from web pages and libraries. Objective 5 will be realized by studying the handouts in the course of RAMS Assessment and Optimisation at NTNU. Besides, supplemented articles and useful insight from professor Anne Barros has contributed valuable inputs in the identification and analysis of the problems in the thesis.

1.5 Structure of the Report

The outline of this report is as follows. Chapter 2 presents a brief survey of the review papers related to prognostic with RUL. Modeling approaches like regression, Brownian motion, and Gamma process are briefly reviewed, advantages and limitations are stated in each section.

(15)

CHAPTER 1. INTRODUCTION 6 Chapter 3 and 4 presents the most important statistical properties of two stochastic process:

Brownian motion and Gamma process. Parameter estimation processes are implemented based on testing dataset and RUL estimation in degradation models are also presented. Chapter 5 studies different monitoring methods and maintenance policies, and discusses how stochastic process can be used to optimize condition-based maintenance. Chapter 6 summarizes what has been done in this thesis and recommendations of future work.

(16)

Chapter 2

Review of the Problem

Maintenance decision making is a very critical step in condition-based maintenance. In order to support maintenance decision in a CBM program, two kinds of techniques can be involved: di- agnostics and prognostics. AsJardine et al.(2006) states, diagnostics focuses on detection, isola- tion and identification of faults when they occur, while prognostic, on the other hand, attempts to predict faults or failures before they occur. From the definitions, we can see that prognostic is a more proactive approach as we can grossly predict and prevent upcoming failures and be prepare for the problems, and thus save extra unplanned maintenance cost. Generally, deteri- oration is monitored during working life of the subject so as to anticipate any faults in time. In the following sections, the writer will review articles that related to degradation modeling and prognostic in CBM decision support.

2.1 Literature Review on Degradation Modeling and Prognostic Approaches

With the rapid advance in prognostic, several review papers specifically on degradation model- ing and maintenance issues have appeared, see for example (Jardine et al.(2006),Si et al.(2011), Ahmadzadeh and Lundberg (2014)). Among these papers, prognostic models can be grossly divided into four categories: physical, experimental, data-driven and hybrid. Most papers on machine prognostics discuss only data driven and model based approaches, and few address

7

(17)

CHAPTER 2. REVIEW OF THE PROBLEM 8 the experience based approach. In this section, the review will be organized according to these categories.

Physical models assume that an accurate mathematical model can be constructed from first principles. They build technically comprehensive theoretical models to describe the physics of the system and failure modes, like crack propagation, wear corrosion and spall growth (Ah- madzadeh and Lundberg(2014)). Li et al.(1999) relates rolling element bearing defect growth rate to the instantaneous defect area size and material constants. They proposes a remaining life adaptation methodology based on mechanistic modeling and parameter tuning.Li and Lee (2005) presents a model-based method to forecast the RUL based on the estimated crack size and dynamic load.Watson et al.(2005) presents RUL prediction of highly dynamic, high power dry clutch system by combining physics-based simulation and wear prediction models. In these paper, the physical models must be configured specifically for the systems being monitored in order to give accurate prediction. The limitations are their high costs and component specialty.

In other words, they cannot be directly applied to other types of components. Additionally, it could be difficult or even impossible to catch some system’s behavior by a mathematical model.

Experimental based methodology utilizes experiments to improve the understanding of the lifetime of components. Given the data and knowledge accumulated from experience, this ap- proach applies probabilistic or stochastic models of degradation of an subject. Lu and Meeker (1993) use fatigue-crack-growth data to develop statistical methods and use Monte Carlo sim- ulation to obtain point estimates and confidence intervals for reliability assessment. Sutrisno et al.(2012) selects an experimental data set from 17 ball bearings from FEMTO-ST Institute.

They develops prognostic algorithms to estimate RUL of the test bearings. Similarly, Le Son et al.(2013) presents a probabilistic method for prognostic applied to the 2008 PHM Conference Challenge data. A stochastic process is proposed to model the deterioration of the components and to estimate the RUL on a case study. Normally, experimental test are designed and set up in order to simulate real working conditions for studied subjects. Since a critical events could be catastrophic, we need to perform further research and combine this method with other ap- proaches to secure the reliability of subjects.

Data-driven approaches use real data (like online gathered data with sensors or operator measures). They approximate and track features which reveal the degradation of components

(18)

CHAPTER 2. REVIEW OF THE PROBLEM 9 and then forecast the global behavior of a system. They can transform high-dimensional noisy data into lower dimensional logical information, and thus easier to be applied in practice.Chin- nam and Baruah(2003) presents a novel method for employing Hidden Markov Model for au- tonomous diagnostics as well as prognostics, as it exploits competitive learning to achieve HMM- based clustering. Gebraeel and Lawley(2008) focuses on the development of a neural network- based degradation model that utilizes condition-based sensory signals to compute and con- tinuously update residual life distribution of partially degraded components. Si et al.(2013) presents a degradation path-dependent approach for RUL estimation through the combination of Bayesian updating and expectation maximization algorithm. However, the drawbacks of data driven approaches are that they highly depend on the quantity and quality of system operational data.

Hybrid approaches combine two or more prognostic methods for the formulation of degra- dation models and data analysis. Such combination could help to extract degradation data, greatly reduce the complexity of calculation or improve precision in predicting RUL. For in- stance,Banjevic and Jardine(2006) calculates the RUL by considering the hazard rate function and Markov property as a stochastic covariate process. Illustration of the main concepts is given using field data from a transmission’s oil analysis histories.Mazhar et al.(2007) proposes a com- prehensive two-step approach for RUL estimation of used components in consumer products.

In the first stage, Weibull analysis is applied to the time-to-failure data to assess the mean life of components. In the second stage, the degradation and condition monitoring data are analyzed by developing an artificial network model. Satish and Sarma(2005) combines neural networks and fuzzy logic and forming a fuzzy back propagation network for identifying the present con- dition of the bearing and estimate the RUL of the motor. Yan and Lee(2007) presents a hybrid method for online assessment and performance prediction of RUL in drilling operations based on vibration signals. Logistic regression analysis combined with maximum likelihood technique is employed to evaluate tool wear condition based on features extracted from vibration signals using wavelet packet decomposition technique. Auto regressive moving average model is then applied to predict RUL.

(19)

CHAPTER 2. REVIEW OF THE PROBLEM 10

2.2 Degradation Modeling

In Fig.1.1in Chapter 1, the writer have already introduced the framework of degradation mod- eling based on RUL estimation. As discussed before, the topic of this thesis would be narrowed down to the first category of heterogeneity (unit-to-unit variability). Variability in the inner structures of the considered systems and the diversity in their working conditions would both contribute to unit-to-unit variability in performed degradation. So it is very easy to think of in- cluding random variables to capture such variability, and leaving the rest of parameters as con- stants describing the homogeneity in system degradation from the same type when modeling degradation process and estimate RUL. So, in this section, the writer would briefly review both two subcategories - random coefficients regression model and stochastic process with random effects. The writer only considers Brownian motion and Gamma processes in stochastic process sections as a study preparation for a further investigation in the following chapters.

2.2.1 Random coefficient regression Models

Lu and Meeker(1993) first described a general nonlinear regression model to characterize the degradation of a population of units as

Y(t)=D(t;φ,θ)+ε(t) (2.1)

whereD(t;φ,θ) is the actual degradation path at timet,φis the fixed regression coefficients for all units,θ is the random effect for individual unit, andε(t) is the random noise described by N(0,σε). Notice that θandεare usually assumed to be independent of each other. With this model, the RUL can be defined as

RULti ={hti:D(ti+hti;φ,θ)≥L|D(ti;φ,θ)<L} (2.2)

whereLis a predefined threshold,htiis the residual useful lifetime given the current component status at time ti. Notice that the underlying assumptions of this model include (Wang et al.

(2000)):

1. the condition of the component deteriorates with operating time and the level of degra-

(20)

CHAPTER 2. REVIEW OF THE PROBLEM 11

dation can be observed at any time;

2. the device being monitored comes from a population of devices, each of which exhibits the same degradation form;

3. the distribution of the random variable across the population of components is known.

Based on these assumptions, researchers have extended and developed the model presented byLu and Meeker(1993) by incorporating some differences. For example,Gebraeel et al.(2005) develop a Bayesian updating methods to construct a closed-form residual-life distribution for a monitored device using sensor based monitoring signals. The focus on only single operating device is quite different from methods of computing RUL for a population of components. Some papers also replace Brownian motion for normally distributed noise variables in the degradation model, which enhance the dynamics of the original regression models.

The random deterioration path model is a very simple model to study, and it is directly re- lated to statistical analysis of deterioration data. However, there are several limitations with this model. First, The fundamental assumption of random deterioration path models about the sample space and sample function of the deterioration process is restrictive when the patterns of some sample deterioration paths are inconsistent with the others due to slight or intensive variations in the environment that an individual asset operates. Also, the assumptions of the the common use of an independent and identical normal noise term are quite restrictive in capturing the temporal uncertainty. Another problem is that one single inspection is sufficient to determine the sample path in the linear random variable degradation and thus having no first hitting time (FHT) motivation. So it cannot model the temporal variability in RUL estimation as argued byPandey and Yuan(2006).

2.2.2 Brownian Motion

Continuous-time Markov process with independent increments such as Brownian motion with drift is also termed as the Gaussian or Wiener process. It is a continuous-time stochastic process with drift parameterµand variance parameterσ. The Brownian motion model has an additive

(21)

CHAPTER 2. REVIEW OF THE PROBLEM 12 effect on degradation process and can be expressed as follows,

X(t)=X(0)+µt+σB(t) (2.3)

whereX(0) is the initial degradation value,µtis the trend andσ2is constant diffusion parame- ter. With this model, the RUL can be defined as

RULti =inf{hti :Y(ti+hti)≥L|Y(ti)<L} (2.4)

The PDF of the FHT of the Brownian motion is the inverse Gaussian distribution. So explicit equations for calculation is available (more discussion in Chapter 3).

Improvements have been made for the application of Brownian motion model in both re- liability engineering and biostatistics in relation to estimating the lifetime. For instance,Tseng and Peng(2004) remedies the weakness of Markov property in Wiener process and proposes an integrated Wiener process to model the cumulative degradation path of the product’s quality characteristic. Peng and Tseng(2009) incorporates the random effect of a drift coefficient and measurement errors into a Wiener process-based degradation process for lifetime analysis.

However, the limitation of Brownian motion deserves a few comments. First, the property of either increasing or decreasing in the context of reliability is quite restrictive, and thus inad- equate in modeling degradation which is monotone. Second, a Brownian motion is a time ho- mogeneous process, but it may not apply to all degradation processes, like fatigue crack which the length grows faster during the course of crack propagation and thus producing time hetero- geneity. Additionally, the variance of the noise term in the Brownian motion is proportional to the time interval, which is a strong requirement few state processes can satisfied.

2.2.3 Gamma Process

Since degradation processes are generally monotonic, it is a nature choice to model them by Gamma process, where degradation occurs gradually over time in a sequence of tiny positive increments. With this model, RUL can be can be defined in a similar way as in the previous

(22)

CHAPTER 2. REVIEW OF THE PROBLEM 13 sections,

RULti={hti :Y(ti+hti)≥L|Y(ti)<L} (2.5) Notice that, the sum of Gamma distributed increments is still a Gamma process. Obviously, the calculation is straightforward and can be figured out using the properties of Gamma pro- cess (see more discussion in Chapter 4). It is very useful in providing optimal maintenance de- cisions. Some extensions were for exampleLawless and Crowder(2004) constructs a tractable gamma process model incorporating a random effect to characterize different degradation rates among individual components.Kuniewski et al.(2009) presents a method to model degradation process which is a combination of two stochastic processes, namely the process of defect initi- ation (non-homogeneous Poisson process) and the process of defect growth (Gamma process).

In a word, Gamma process is chosen for its relatively straightforward calculation and ability to capture temporal variability.

But it should be remembered that Gamma process appears only appropriate to model degra- dation by a strictly monotonic process. Also, it is a challenge to reasonably choose random parameters and their distributions which both capture the unit-to-unit variability and simplify RUL calculation. Given these reasons, some modified Gamma processes like semi-parametric gamma process should be further studied to improve their modeling ability.

(23)

Chapter 3

Brownian Motion Model

Brownian motion is a continuous stochastic process that is widely used for modeling random behavior evolves over time. Today, the process and its many generalizations and extensions occur in diverse areas of both pure and applied science, such as mathematical statistics, man- agement science,economics, communication theory and biology etc.

The purpose of this section is to review the some mathematical properties and understand the advantages of Brownian motion. The merits could possibly be extended to the understand- ing and application of other Markov stochastic processes, such as Gamma process. This chapter would first present some mathematical properties of standard Brownian motion (i.e. Wiener process), and then move to Brownian motion with drift. Besides, methods for Brownian motion path simulation, parameter estimation and RUL distribution estimation are presented in the following sections.

3.1 Standard Brownian Motion - Wiener Process

3.1.1 Introduction and Definition

Consider a symmetric random walk, which in each time unit is equally likely to take a unit step either up or down . If this process is speed up by taking smaller and smaller steps in smaller and smaller time intervals, we would obtain Brownian motion once we go the limit in the right manner.

14

(24)

CHAPTER 3. BROWNIAN MOTION MODEL 15 For ease, we start with the assumptionX(0)=0. We also assume that there is no "drift" to the process, i.e. E(X(t))=0. Let {ξn}n∈Nbe a sequence of independent, identically distributed (i.i.d) random variables such that

E(ξn)=0,

Var(ξn)=(ξ2n)=1, (3.1)

Notice thatξnsatisfies

ξn=

+1, if theith step of lengthξnis up

−1, if it is down

(3.2)

ThenX(n) is the instantaneous position of a random walk on the integersZ. It can be viewed as a function of the discrete timenas follows

X(0)=0, X(n)=

n

X

i=1

ξi, (3.3)

We will then rescale both time and space in order to construct a random function on t ∈ [0,+∞) and taking values inR. Recall that the Central Limit Theorem asserts that

X(N)

pNN(0, 1) (3.4)

in distribution as N → ∞. This suggest to rescale X(n) and construct a constant random functionXN(t) ont∈[0,+∞) by letting

XN(t)=SbN tc

pN (3.5)

wherebN tcrepresents the largest integer less thanN t.

The function is illustrated in Fig.3.1. for different resolutionN. From the figure, we can see that when N increases, the discrete stair process will approach a continuous-time and space process. By the Central Limit Theorem, the distribution of XN(1)= pSNN approaches a normal

(25)

CHAPTER 3. BROWNIAN MOTION MODEL 16

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 100-step wiener process

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 1000-step wiener process

Figure 3.1: Realizations ofXN(t) (i.e.E(X(t))=0) forN=100(Left),N=1000(Right).

distribution with mean 0 and variance 1. Similarly, the distribution ofXN(t) approaches a nor- mal distribution with mean 0 and variance t. It can be shown that asN→ ∞,XN(t) converges in distribution to Brownian motion, that isXN dX, whereX(·) is the Brownian motion.ddenotes that the random process on both sides of the equality have the same distribution. The proof of such theorem is not presented here and we just study its properties in the following sections.

Recall Eq.3.4, we can construct Eq.3.5like this

X(t)=XN(t)=SbN tc pN

=SbN tc pN

pbN tc pn

d p

t N(0, 1)=d N(0,t) (3.6)

Besides, it is easy to prove that increments ofX(t) are normally distributed if we take each

t time unit from the above equation. Standard Brownian motion is a Brownian motion with µ=0,σ2=1. We can also speak of a general Brownian motion withµ6=0,σ26=1 based on the same idea. The expressions can be referred to Summary.

In a word, the Brownian motionX(t) is essentially the accumulation of a series of normally distributed random variables. Between each time interval, the increment of the motion are nor-

(26)

CHAPTER 3. BROWNIAN MOTION MODEL 17 mally distributed with mean equals 0 and variance equals σ2t. As time goes, the variances of these normally distributed random variables increase meaning that it is more uncertain to predict the value of the process after a longer period.

Summary

For∆X =X(t+∆t)X(t)∼σp

t N(0, 1), we have

E(∆X) = 0

Var(∆X) = σ2∆t (3.7)

E(X(t)) = E(X(t)−X(0))=0

Var(X(t)) = Var(X(t)−X(0))=n·σ2t=σ2t

(3.8)

Notice here the variance is additive because any pair of∆X(t)i are assumed to be indepen- dent. Asn→ ∞,∆tconverges to 0 and is normally denoted asd t, which means an infinitesimal time interval. Correspondingly,∆X is re-denoted asd X. However, it is not very precise in math- ematical sense to denoted tandd Xin this way here. We mention it as a notation just to simplify the process of turning discrete stochastic process to a series of continuous normal distributions, i.e., a Brownian motion.

3.1.2 Properties of Standard Brownian Motion

When we speak of a Brownian motion withσ2=1, andX(0)=0, it turns to be a standard Brown- ian motion; otherwise, it is a general Brownian motion. Therefore, a standard Brownian motion X(t) is a stochastic process with the following properties (Lawler(2006)).

1. Initial conditionX(0)=0, the process start at 0;

2. Independence of incrementsFor anys1t1s2t2≤ · · · ≤sntn, the random variables Xt1Xs1, . . . ,XtnXsn are independent;

3. Normal incrementsFor anys<t, the random variableXtXshas a normal distribution with mean 0 and variancets. If we makes=0, thenX(t)−X(0) hasN(0,t) distribution;

(27)

CHAPTER 3. BROWNIAN MOTION MODEL 18 4. Continuity of pathsX(t),t≥0 are continuous functions oft.

3.1.3 Properties of Standard Brownian Path

Given stochastic differential equations (SDEs)

d y(t)=a[y(t),t]d t+b[y(t),t]d x(t) (3.9)

a general representation of Brownian sample path can be constructed as follows,

d y(t)=µd t+σd x(t) (3.10)

whered x(t) is a standard Brownian motion (i.e. Wiener process). With initial condition X(0)=0, scale parameterµ=0 and shape parameterσ=1, Eq.3.10describes the evolution of a standard Brownian motion with meanµt=0 and varianceσ2t=t. So the standard Brownian motion path can also be written as

d y(t)=d x(t) (3.11)

Y(t)=X(t)=X(0)+µt+σ(X(t)−X(0)) (3.12)

The density function ofY(t) given thatY(0)=0 is the transition probability density function of Brownian motion. Letpt(x,y) denote the transition densities fromxtoy,

pt(x,y)= 1

p2πte−(y−x)22t (3.13)

According toKlebaner et al.(2005),X(t) as functions oft have the following properties. Al- most every sample pathX(t), 0≤tT

1. is a continuous function oft;

2. is not monotone in any interval, no matter how small the interval is;

3. is not differentiable at any point;

(28)

CHAPTER 3. BROWNIAN MOTION MODEL 19

4. has infinite variation on any interval, no matter how small it is;

5. has quadratic variation on [0,t] equal tot, for anyt.

The properties are not proved here because they can be found in many books on stochastic processes, likeLawler(2006) andKlebaner et al.(2005).

Summary of Standard Brownian Motion For∆X =X(t+∆t)X(t)∼p

t N(0, 1), we have

E(∆X) = 0

Var(∆X) = ∆t (3.14)

E(X(t)) = E(X(t)−X(0))=0

Var(X(t)) = Var(X(t)−X(0))=n·∆t=t

(3.15)

3.2 Brownian Motion with Linear Drift

3.2.1 Properties of Brownian Motion with Linear Drift

Similar to what have been presented in Section3.1.2, Brownian motion with Linear driftY(t) possesses the following properties

1. Initial conditionY(0)=x, the process start atx;

2. Independence of incrementsFor anys1t1s2t2≤ · · · ≤sntn, the random variables Yt1Ys1, . . . ,YtnYsnare independent;

3. Normal incrementsFor anys<t, the random variableYtYshas a normal distribution with mean (t−s)µand variance (t−s)σ2. If we makes=0, thenX(t)−X(0) hasN(µt,σ2t) distribution;

4. Continuity of pathsY(t),t≥0 are continuous functions oft.

(29)

CHAPTER 3. BROWNIAN MOTION MODEL 20

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 1000-step wiener process and its mean

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 1000-step wiener process and its mean

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 1000-step wiener process and its mean

Mu 0, Sigma 3

0 0.2 0.4 0.6 0.8 1

-6 -4 -2 0 2 4

6 1000-step wiener process and its mean

Figure 3.2: Four realizations of standard Brownian paths.

Refer to Eq.3.10, the Brownian motion path with drift is

d y(t)=µd t+σd x(t) (3.16)

where initial conditionY(0)=x, and it describes the path of a Brownian motion with mean µtand varianceσ2t. Therefore, it can also be expresses as

Y(t)=x+µt+X(t) (3.17)

whereX(t) is a zero drift Brownian motion with varianceσ2starting at 0.

Notice that the motionY(t) consists of a linear motion in the directionµwith random fluc- tuations related toσ.

(30)

CHAPTER 3. BROWNIAN MOTION MODEL 21

Mu 0.1, Sigma 1

0 0.2 0.4 0.6 0.8 1

-1.5 -1 -0.5 0 0.5 1 1.5 2

1000-step wiener process

Mu 1, Sigma 0.1

0 0.2 0.4 0.6 0.8 1

0 0.2 0.4 0.6 0.8 1

1.2 1000-step wiener process

Figure 3.3: Realizations of differentµ,σin Brownian paths.

The unconditional probability density function at a fixed timetgivenY(0)=xis

pt(x,y)= 1 p2πσ2te−|

yx−µt|2

2σ2t (3.18)

Summary of Brownian Motion with Linear Drift For∆Y =Y(t+∆t)Y(t)∼p

t N(µ,σ2), we have

E(∆Y) = µ∆t Var(∆Y) = σ2t

(3.19)

E(Y(t)) = E(Y(t)−Y(0))=x+µt

Var(Y(t)) = Var(Y(t)−Y(0))=n·σ2t=σ2t

(3.20)

3.2.2 Some Graphs and Discussions

Graphs of some examples of Brownian motion are given in order to visualize the processes and the equations presented in the Summary.

Example 1. Randomness Suppose we simulate a Brownian motion process by Matlab with meanµ=0 and varianceσ=0.1. The overall time is 1 and there are 1000 steps. Initial condition

(31)

CHAPTER 3. BROWNIAN MOTION MODEL 22 isY(0)=x=0. Four realizations of such Brownian motion are shown in Fig.3.2. We setE(Y(t))= 0 which represents a process with no drift (i.e. zero mean) and governed by pure randomness (White Noise). However, it appears that the paths has regions where motions looks like they has

"trends".

Example 2. Drift and Noise Suppose we simulate two Brownian motion paths with initial conditionY(0)=x =0. The overall time scale is 1 and there are 1000 steps. We suppose the first path is with parameterµ=0.1,σ=1 (See Fig.3.3, left), while the second path is withµ= 1,σ=0.1(See Fig.3.3, right). These figures demonstrate that when scale parameterµis small in comparison with shape parameter σ, drift has a greater impact on the Brownian process;

if shape parameterσis small in comparison with scale parameterµ, then noise dominates in the behaviour of the Brownian process. Recall from Eq. 3.17, X(t) can be considered as "Pure Noise", while (x+µt) is the "Actual Degradation Path" of interest.

Example 3. Effect ofµandσ Suppose we simulate four Brownian motion paths with initial conditionY(0)=x=0. The overall time scale is 1 and there are 1000 steps. We first keep the value ofσ=1 in all four realizations and only changeµto see its impact on process. The paths are with parametersµ=1, 2, 3, 4. Fig. 3.4shows 200 realizations and their means (blue lines) in each graph. Recall Eq.3.20and we can see from the figures that the slopeE(Y(t)) increases when µincreases. It is not clear from the figure whether the processes are more concentrated, spread out or just keep the same, but from Eq.3.20, the variance would not be affected by changes inµ and keep all the same.

Then, we keep the value ofµ=1 in all four realizations and only changeσto see its impact on process. The paths are with parameters µ=1, 2, 3, 4. Fig. 3.5shows 200 realizations and their means(blue lines) in each graph. Recall Eq.3.20and we can see from the figures that the slope (E(Y(t))) keeps unchanged whenµincreases. But it is very obvious from the figure that the processes are more spread out. It demonstrates Eq. 3.20that the variance of process increases whenσincreases.

(32)

CHAPTER 3. BROWNIAN MOTION MODEL 23

Figure 3.4: Four groups of realizations with increasingµ.

3.2.3 First Hitting Time and RUL Distribution

One of the purposes of studying degradation model is to estimate the residual useful lifetime and its distribution function. If previous data before a monitoring time t are available, how much time would the component still have to be useful? The first hitting time is one of the challenges when we transfer this situation into Brownian motion models. From the previous definitions and figures, we can easily find out that Brownian motion is a non-monotone process.

Given a levelL, the time of hitting this level could be more than once due to the random nature of Brownian motion process. The following contents provide the general computation of a first passage time hits levelL, hits timeTLand RUL distribution.

Hitting LevelL SupposeTLis the first time Brownian motionY(t) hits levelL. Let

(33)

CHAPTER 3. BROWNIAN MOTION MODEL 24

Figure 3.5: Four groups of realizations with increasingσ.

M(t)=max

0≤s≤tY(s) (3.21)

TL=inf{t>0 :Y(t)≥L} (3.22)

If the maximum at timet is greater thanL, then the Brownian motion took valueLat some time beforet. Meanwhile, if Brownian motion took valueLat a certain time beforet, then the maximum value will be at leastL. In mathematical form, the above argument could be written as follows,

{Y(t)≥L}⊂{TLt} (3.23)

(34)

CHAPTER 3. BROWNIAN MOTION MODEL 25 Therefore, it is nature to say that

P(Y(t)≥L)=P(Y(t)≥L,TLt) (3.24)

SinceY(TL)=L,

P(Y(t)≥L)=P(TLt,Y(TL+(t−TL))−Y(TL)≥0) (3.25)

By TheoremB.1in appendix,TLis a finite stopping time, and by the strong Markov property (B.2), the random variable ˆY(s)=Y(TL+s)−Y(TL) is independent ofFTL and has a Normal distribution, so

P(Y(t)≥L)=P(TLt, ˆY(t−TL)≥0) (3.26) It is easy to show thatTLt andM(t)L are two equal events. Ifs is independent ofTL, then

P(TLt, ˆY(s)≥0)=P(TLt)P( ˆY(s≥0))

=P(TLt)1 2

=P(M(t)≥L)1

2 (3.27)

That is to say for anyL>0,

Px(M(t)≥L)=2P0(Y(t)≥L)

=2 µ

1−Φ

µLxµt σp

t

¶¶

(3.28)

wherePxdenotes the probability of events when process starts atx,Φ(·) stands for the standard Normal distribution function. The last part of equation can be explained by central limit theo- rem. But in Eq.3.26s=tTL, and is clearly dependent onTL. It is not easy to transfer Eq.3.26to Eq.3.28. Detail proof is not discussed here, but can be found in other books.

(35)

CHAPTER 3. BROWNIAN MOTION MODEL 26 First Passage TimeTL To derive the distribution of first passage timeTL, we should first mention inverse Gaussian distribution (IG). It is a two-parameter continuous distribution given by its density function

f(x;µ,λ)=

· λ 2πx3

¸1/2

exp−λ(x−µ)2

2x (3.29)

forx>0, whereµ>0 is the mean andλ>0 is the shape parameter. For a random variable X with inverse Gaussian distribution we writeXIG(µ,λ).

According toPärna (2014), the inverse Gaussian distribution describes the distribution of the time a Brownian motion with positive drift takes to reach a given positive level. As what has been denoted before, letTL be the first passage time for a fixed levelL>0 byYt. ThenTL has inverse Gaussian distribution, TLIG(Lµ,σL22). So put these parameters into Eq.3.29, we have probability density function

f(x;µ,σ)= L

p2πµ3x2exp µ

−(L−µx)2 2σ2x

(3.30)

Therefore, the first passage timeTLsatisfies the following function

F(x;µ,σ)=P(TLt)= Z t

0

L

p2πµ3x2exp µ

−(L−µx)2 2σ2x

d x (3.31)

RUL Distribution Given Eq.3.30, it is easy to obtain residual useful liftime distribution.

Observing the process at timetis at positiony(t), the probability that residual useful lifetime is less than a predefined periodhis

P(RU Ly(t)h)= Z h

0

Ly(t) p2πµ3x2exp

µ

−(L−y(t)−µx)22x

d x (3.32)

3.3 Philosophy of RUL Estimation

The purpose of this section to estimate remaining useful lifetime (RUL) using Brownian motion process to model the degradation. First, assumptions and an example training data set are pre- sented. Then, the parameters in the degradation model are estimated with the training data set.

And finally the RUL distribution is estimated proposed for the units of testing data.

(36)

CHAPTER 3. BROWNIAN MOTION MODEL 27

3.3.1 Assumptions and Training Data Set

In this example, the Brownian motion process is used for degradation modeling correspond- ing to training data set. Assume thatY is an degradation indicator forN =100 independent identical tested items, andYi,j denotes the degradation measurements of theith items at time j, where i =1, 2, . . .N and j =1, 2, . . .r, r is the last observation time. All the degradation re- alizations are based on homogeneous Brownian motion with drift with two parameters (µ,σ), and that are the same for all the items. According to Eq.3.19, each increment of degradation

Yi,j =Yi,j+1Yi,j of each items follows a normal distribution N(µ∆ti,j,σ2ti,j), whereµ= 6,σ=4. There are 1000 observations within the time frameT =1, and the initial condition is Y(0)=x=1. Let the time step equals 1/1000. The degradation paths for allN =100 units and certain units are illustrated in Fig.3.6. The following parameter estimation, process verification and RUL distribution are all based on these built testing dataset.

3.3.2 Parameter Estimation

Given degradation pathYi jbased on Brownian motion process, each increment is∆Yi,j =Yi,j+1Yi,jN(µ∆ti,j,σ2ti,j) is normally distributed for all identically and independent components.

Recall the transition density function of Brownian motion process in Eq.3.18, then

f(µ∆ti,j2∆ti,j)(∆Yi,j)= 1 q

2πσ2ti,j e

(Yi,j−µ∆ti,j)2

2σ2ti,j (3.33)

SinceY is Markovian, the maximum likelihood estimator (MLE) ofθ=(µ,σ) can be calcu- lated once transition density function ofY are known. Consider the degradation measurements for itemi,∆Yi=(∆Yi,1,∆Yi,2, . . . ,∆Yi,r). For itemi, the likelihood function for itemi is

(37)

CHAPTER 3. BROWNIAN MOTION MODEL 28

0 0.2 0.4 0.6 0.8 1

-5 0 5 10 15

20 Degradation paths for all units

0 0.2 0.4 0.6 0.8 1

-5 0 5 10 15

20 Unit 1

0 0.2 0.4 0.6 0.8 1

-5 0 5 10 15

20 Unit 50

Mu 6, Sigma 4

0 0.2 0.4 0.6 0.8 1

-5 0 5 10 15

20 Unit 100

Figure 3.6: Degradation paths of all units and certain units.

Li(θ)=fi(∆Yi)=fi(∆Yi,1,∆Yi,2, . . . ,∆Yi,r|µ,σ)

=

r

Y

j=1

fi(∆Yi,j|µ,σ)

=

r

Y

j=1

1 q

2πσ2ti,j e

(Yi,j−µ∆ti,j)2

2σ2ti,j (3.34)

(38)

CHAPTER 3. BROWNIAN MOTION MODEL 29 Then for theith item, the log-likelihood is given by

li(θ)=lnLi(θ)=ln

 Yr j=1

1 q

2πσ2ti,j e

(Yi,j−µ∆ti,j)2 2σ2ti,j

 (3.35)

Since the measurementsYi,j are independent

l(θ)=ln(∆Y1,∆Y2, . . . ,∆YN)

=

N

X

i=1

ln(fi(∆Yi,1,∆Yi,2, . . . ,∆Yi,r))

=

N

X

i=1

ln

r

Y

j=1

1 q

2πσ2ti,j e

(Yi,j−µ∆ti,j)2 2σ2ti,j

 (3.36)

that is to sayl(θ)=PN

i=1li(θ), wherefi/li(θ) is the probability density function /log-likelihood of increments corresponding to each item, and f/l(θ) is the probability density function /log- likelihood of all increments.

The maximum likelihood estimator ˆθ=( ˆµ, ˆσ) are found by maximizingl(θ). In practice it is done by taking the partial derivative of the log-likelihood function of Eq. 3.36with respect toµ andσ. This gives equations

∂l(θ)

∂µ = XN i=1

Xr j=1

Yi,jµ∆ti,j

σ2 =0 (3.37)

∂l(θ)

∂σ = −r N σ +

N

X

i=1 r

X

j=1

(∆Yi,jµ∆ti,j)2

σ3ti,j =0 (3.38)

So the maximum likelihood estimator forθ=(µ,σ) is

µˆ= PN

i=1

Pr

j=1Yi,j PN

i=1

Pr

j=1∆ti,j

(3.39)

σˆ= v u u t

1 r N

N

X

i=1 r

X

j=1

(∆Yi,jµ∆ti,j)2

ti,j (3.40)

Eq.3.40can be substituted and solved by Eq.3.39.

Referanser

RELATERTE DOKUMENTER