Reliability measures for indexed semi-Markov chains applied to wind energy production
RReliability measures for indexed semi-Markov chainsapplied to wind energy production
Guglielmo D’Amico
Dipartimento di Farmacia, Universit`a ‘G. D’Annunzio’ di Chieti-Pescara, 66013 Chieti,Italy
Filippo Petroni
Dipartimento di Scienze Economiche ed Aziendali, Universit`a degli studi di Cagliari,09123 Cagliari, Italy
Flavio Prattico
Dipartimento di Ingegneria Industriale e dell’Informazione e di Economia, Universit`adegli studi dell’Aquila, 67100 L’Aquila, Italy
Abstract
The computation of the dependability measures is a crucial point in theplanning and development of a wind farm. In this paper we address the issueof energy production by wind turbine by using an indexed semi-Markov chainas a model of wind speed. We present the mathematical model, we describethe data and technical characteristics of a commercial wind turbine (AirconHAWT-10kW). We show how to compute some of the main dependabilitymeasures such as reliability, availability and maintainability functions. Wecompare the results of the model with real energy production obtained fromdata available in the Lastem station (Italy) and sampled every 10 minutes.
Keywords: semi-Markov chains; synthetic time series; dependabilityanalysis
1. Introduction
Wind is one of the most important renewable energy sources. Windenergy is produced by converting the kinetic energy of wind into electricalenergy by means of a generator. For this reason it is important to dispose
Preprint submitted to Probabilistic Engineering Mechanics October 17, 2018 a r X i v : . [ phy s i c s . d a t a - a n ] N ov f an efficient stochastic model for wind speed changes. In a recent paper [1]the authors showed that an indexed semi-Markov chain (ISMC) reproducesalmost exactly the statistical features of wind speed. In particular, it wasshown that density and autocorrelation functions of a time series of realwind speed and those obtained by the ISMC model through Monte Carlosimulation were almost indistinguishable.In this work we use the ISMC model to compute dependability measuresas availability, reliability and maintainability functions for the index semi-Markov process. These indicators give important information on the feasi-bility of the investment in a wind farm by giving the possibility to quantifythe uncertainty in the wind energy production.Another important aspect is related to the location of the wind farm. Infact, today many wind farm are built offshore for different reasons: the windspeed is more powerful and constant due to the absence of obstacles, the vi-sual, environmental and acoustic impact is cut down. The maintenance cost,instead, is higher than the onshore wind farm. A good stochastic model canhelp the planning of preventive maintenance suggesting when it is suitableto do the maintenance operation.The results presented here are new and generalize some of the results ob-tained for semi-Markov chain(see [2, 3, 4, 5, 6, 7, 8]). The model generalizesalso Markov chains and renewal models.We apply our model to a real case of energy production. For this reasonwe choose a commercial wind turbine, a 10 kW Aricon HAWT assumed tobe installed at the station of L.S.I -Lastem which is situated in Italy.The paper is organized as follows. Section 2 presents some definitionsand notation on the indexed semi-Markov chain. Section 3 describes thedatabase and the technical characteristics of the commercial wind turbine.Section 4 shows the way in which it is possible to compute the dependabilitymeasures via kernel transformations and the value computed on the real dataand on the synthetic data are compared. In the last section some concludingremarks and possible extensions are presented.
2. The indexed semi-Markov process
In this section, we give a short review of a particular indexed semi-Markovmodel advanced in [1] as a novel synthetic time series generation method for2ind speed data which is able to capture the persistence in the wind speeddata and we give new probabilistic result about the transition probabilityfunction.Let us consider the stochastic process J − ( m +1) , J − m , J − ( m − , ..., J − , J , J , ... (1)with a finite state space E = { , , ..., S } . In our framework the randomvariable J n describes the wind speed at the n -th transition, that is at the n -th change of speed value and, the state space E , describes the discretizedwind speeds, see Section 3 for a specific choice of E .Let us also consider the stochastic process T − ( m +1) , T − m , T − ( m − , ..., T − , T , T , ... (2)with values in IN. The random variable T n describes the time in which the n -th change of value of the wind speed occurs. We denote the stochasticprocess { X n } n ∈ IN of the sojourn times in wind speed state J n − before the n th jump. Thus we have for all n ∈ IN X n = T n +1 − T n .Let us consider another stochastic process V − ( m +1) , V − m , V − ( m − , ..., V − , V , V , ... (3)with values in IR. The random variable V n describes the value of the indexprocess at the n th transition: V mn = 1 T n − T n − ( m +1) m (cid:88) k =0 X n − − k (cid:88) s =1 f ( J n − − k , s ) , (4)where f : E × IN → IR is a generic function and V m − ( m +1) , ..., V m are knownand non-random.The function f depends on the state of the system J n − − k and on thetime s .The process V mn can be interpreted as a moving average of the accumu-lated reward process with the function f as a measure of the permanencereward per unit time.It should be noted that the order of the moving average is on the numberof transitions m + 1 which is fixed. Anyway, the moving average is executedon time windows of variable length T n − T n − ( m +1) = (cid:80) nk = n − ( m +1) X k because3ach transition has a random sojourn time X k of permanence in state J k − before the next jump.The indexed model is fully specified once the dependence structure be-tween the variables is assumed. Toward this end we adopt the followingassumption: P [ J n +1 = j, T n +1 − T n ≤ t | σ ( J h , T h ) , h = − m, ..., , ..., n, J n = i, V mn = v ]= P [ J n +1 = j, T n +1 − T n ≤ t | J n = i, V mn = v ] =: Q mij ( v ; t ) , (5)where σ ( J h , T h ) , h ≤ n represents the set of past values, of the processes( J, T ). The matrix of functions Q m ( v ; t ) = ( Q mij ( v ; t )) i,j ∈ E is called indexedsemi - M arkov kernel .The joint process ( J n , T n ), which is embedded in the indexed semi-Markovkernel, depends on the moving average process V mn , the latter acts as astochastic index. Moreover, the index process V mn depends on ( J n , T n ) throughthe functional relationship (4).Observe that if P [ J n +1 = j, T n +1 − T n ≤ t | J n = i, V λn = v ] = P [ J n +1 = j, T n +1 − T n ≤ t | J n = i ] , (6)for all values v ∈ IR of the index process, then the indexed semi-Markovkernel degenerates in an ordinary semi-Markov kernel and the model becomesequivalent to classical semi-Markov chain models as presented for examplein [9, 10, 11, 12, 13, 14]. The dependence of the process ( J n , T n ) on the newvariable V mn is introduced in order to capture the effect of the past transitionson the future ones for those processes which are strongly autocorrelated.One of the main problem is the proposal of a particular choice of theindex process which is useful in describing the wind speed process. To thisend we need to choose a specific form of the function f .The choice is motivated by some physical reasons and by model simplicity.Let us briefly remind that wind speed data are long range positivelyautocorrelated. This implies that there are periods of high and low speed.Motivated by the empirical facts in [1] we supposed that also the transitionprobabilities from current wind speed J n to the next one J n +1 depends onwhether the wind is, on average, in a high speed period or in a low one. Wethen fixed the function f to be the wind speed itself, i.e. f ( J n − − k , s ) = J n − − k , (7)4or all s ∈ IN. Consequently, substituting (7) in equation (4) and consideringthat J n − − k is constant in s we obtain V mn = 1 T n − T n − ( m +1) m (cid:88) k =0 J n − − k · X n − − k = m (cid:88) k =0 J n − − k · (cid:18) T n − k − T n − − k T n − T n − ( m +1) (cid:19) . (8)In this simple case the index process ( V mn ) expresses a moving averageof order m + 1 executed on the series of the wind speed values ( J n − − k )with weights given by the fractions of sojourn times in that wind speed( T n − k − T n − − k ), with respect to the interval time on which the average isexecuted ( T n − T n − ( m +1) ).The consequence of the choice (7) is that the assumption (5) now statesthat: Q mij ( v ; t ) = P [ J n +1 = j, T n +1 − T n ≤ t | J n = i, V mn = m (cid:88) k =0 J n − − k (cid:18) T n − k − T n − − k T n − T n − ( m +1) (cid:19) = v ] . (9)Relation (9) asserts that the knowledge of the last wind speed value ( J n = i ) and of the weighted moving average ( V mn = v ) of order m + 1 of past windspeeds suffices to give the conditional distribution of the couple J n +1 , T n +1 − T n whatever the values of the past variables might be. Essentially we considerthat the average of the past m speeds contains most of the information neededto establish the probability of the next wind speed transition. We will showlater that, with this assumption, the model is able to capture the temporaldependence structure of real data.We introduce now auxiliary probabilities which are helpful in the sequelof the analysis. Denote by p mij ( v ) := P [ J n +1 = j | J n = i, V mn = v ] . They represent the transition probabilities of the embedded indexed Markovchain. More precisely p mij ( v ) denotes the probability that the next transitionis into wind speed j given that at current time the wind speed process enteredstate i and the index process had value v . It is simple to realize that p mij ( v ) = lim t →∞ Q mij ( v ; t ) . (10)Let H mi ( v ; · ) be the sojourn time cumulative distribution in wind speedstate i ∈ E : H mi ( v ; t ) := P [ T n +1 − T n ≤ t | J n = i, V mn = v ] = (cid:88) j ∈ E Q mij ( v ; t ) . (11)5t expresses the probability to change the actual wind speed i in a timeless or equal to t given the indexed process has value v .It is useful to consider also the conditional waiting time distribution func-tion G which expresses the following probability: G mij ( v ; t ) := P [ T n +1 − T n ≤ t | J n = i, J n +1 = j, V mn = v ] . (12)It is simple to establish that G mij ( v ; t ) = (cid:40) Q mij ( v ; t ) p mij ( v ) if p mij ( v ) (cid:54) = 01 if p mij ( v ) = 0 . (13)To describe the behavior of our model at whatever time t ∈ IN we needto define additional stochastic processes.Given the three-dimensional process { J n , T n , V mn } and the indexed semi-Markov kernel Q m ( v ; t ), we define by N ( t ) = sup { n ∈ N : T n ≤ t } ; Z ( t ) = J N ( t ) ; V m ( t ) = 1 t − T ( N ( t ) − θ ) − m m (cid:88) k =0 J ( N ( t ) − θ ) − k · ( t ∧ T ( N ( t ) − θ )+1 − k − T ( N ( t ) − θ ) − k )(14)where T N ( t ) ≤ t < T N ( t )+1 and θ = 1 { t = T N ( t ) } .The stochastic processes defined in (14) represent the number of transi-tions up to time t , the state of the system (wind speed) at time t and thevalue of the index process (moving average of function of wind speed) up to t , respectively. We refer to Z ( t ) as an indexed semi-Markov chain.The process V m ( t ) extends the process V mn because the time t can be atransition or a waiting time. It is simple to realize that, ∀ m , if t = T n wehave that V m ( t ) = V mn .Let us introduce the stochastic process B ( t ) = t − T N ( t ) , (15)which is called backward recurrence time process and denotes the time sincethe last transition.This process is very important in a semi-Markovian framework and it iswell known that the transition probabilities of a semi-Markov process depend6n the value of the recurrence time process, see e.g. [6, 15].This is due to the fact that the conditional waiting time distributionfunctions (12) can be of any type and then, also no memoryless distributionscan be used. In this case, the time since the last transition of the system(backward value) influences the system’s transition probabilities; this is theso called duration effect.To properly assess the probabilistic behavior of the system, we need tointroduce the transition probability function. To this end, it is useful tointroduce the following notation: J − m − = ( J , J − , ..., J − m − ) , T − m − = ( T , T − , ..., T − m − ) . With the equality J − m − = j − m − we denote the fact that( J = j , J − = j , ..., J − m − = j − m − ) , and similarly with T − m − = t − m − we indicate that( T = t , T − = t , ..., T − m − = t − m − ) . Finally by ( J , T ) − m − we denote the couple of ordered vectors ( J − m − , T − m − ). Definition 1.
The transition probability function of the indexed semi-Markovchain Z ( t ) is the function b φ ( i − m − ; j ) ( t − m − ; u, t ) defined by b φ ( i − m − ; j ) ( t − m − ; u, t ) := b φ ( i − m − ,i − m ,...i ; j ) ( t − m − , t − m , ..., t ; u, t ) , (16) where, i , i − , ..., i − m − ∈ E and t , t − , ..., t − m − ∈ IN . The transition probability function (16) expresses the probability the ISMCoccupies state j at time t given that at current time t + u it is in state i where it entered with last transition at time t having previously visited thestates i − − m − at times t − − m − .We call the transition probabilities (16) transition probabilities with ini-tial backward times for the ISMC model. In fact, at the current time t + u the backward recurrence time process assumes value: B ( t + u ) = t + u − T N ( t + u ) = t + u − T = t + u − t = u.
7f we set u = 0 in (16) we have the probability φ ( i − m − ; j ) ( t − m − ; t ) := P [ Z ( t ) = j | J = i , ..., J − m − = i − m − , T = t , ..., T − m − = t − m − ] . (17)which denotes the probability that the ISMC occupies state j at time t given that at current time t it is entered into state i having previouslyvisited the states i − − m − at times t − − m − .The following result consists in a recursive formula for computing thetransition function b φ ( i − m − ; j ) ( t − m − ; u, t ) of the ISMC Z ( t ). Proposition 2.
The probabilities b φ ( i − m − ; j ) ( t − m − ; u, t ) verify the followingequation: b φ ( i − m − ; j ) ( t − m − ; u, t )= δ i j (cid:32) − H mi (cid:16) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; t − t (cid:17)(cid:33)(cid:32) − H mi (cid:16) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; u (cid:17)(cid:33) + (cid:88) i ∈ E t (cid:88) t = t + u +1 q i i (cid:0) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; t − t (cid:1)(cid:32) − H mi (cid:16) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; u (cid:17)(cid:33) b φ ( i − m ; j ) ( t − m ; t − t ) , (18) where δ i j is the Kronecker delta and q ir ( v ; s ) = Q ir ( v ; s ) − Q ir ( v ; s − . Proof
First of all let us compute the value of the index process V m ( t )given the information set { ( J , T ) − m − = ( i − m − ; t − m − ) } . Because t = T is a transition time, we have that θ = 1. Moreover T − m − = t − m − and ∀ kT − k = t − k < t . Then, we have that V m ( t ) = m (cid:88) k =0 m (cid:88) k =0 i − k − t − k − t − k − t − t − m − . (19)Now, being the events { T > t } and { T ≤ t } disjoint, it follows that P [ Z ( t ) = j | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= P [ Z ( t ) = j, T > t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]+ P [ Z ( t ) = j, T ≤ t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ] . (20)8bserve that P [ Z ( t ) = j, T > t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= P [ T > t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ] · P [ Z ( t ) = j | T > t, T > t + u, ( J , T ) − m − = ( i , t ) − m − ] . (21)The first factor on the right hand side of (21) is P [ T > t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= 1 − H mi (cid:16) (cid:80) mk =0 i − k − t − k − t − k − t − t − m − ; t − t (cid:17) − H mi (cid:16) (cid:80) mk =0 i − k − t − k − t − k − t − t − m − ; u (cid:17) . (22)The second factor on the right hand side of (21) is simply P [ Z ( t ) = j | T > t, T > t + u, ( J , T ) − m − = ( i , t ) − m − ] = δ i j , (23)because being T > t the time of next transition exceeds t and, therefore, upto t the process remains in state i . Consequently the probability is one ifand only if j = i otherwise it will be equal to zero.Now let us consider the computation of the second addend on the righthand side of formula (20): P [ Z ( t ) = j, T ≤ t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= (cid:88) i ∈ E t (cid:88) t = t + u +1 P [ Z ( t ) = j, J = i , T = t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= (cid:88) i ∈ E t (cid:88) t = t + u +1 P [ Z ( t ) = j | J = i , T = t , T > t + u, ( J , T ) − m − = ( i , t ) − m − ] × P [ J = i , T = t | T > t + u, ( J , T ) − m − = ( i , t ) − m − ]= (cid:88) i ∈ E t (cid:88) t = t + u +1 b φ ( i − m ; j ) ( t − m ; t − t ) × P [ J = i , T = t | ( J , T ) − m − = ( i , t ) − m − ] P [ T > t + u | ( J , T ) − m − = ( i , t ) − m − ] (24)9 (cid:88) i ∈ E t (cid:88) t = t + u +1 q i i (cid:0) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; t − t (cid:1)(cid:32) − H mi (cid:16) (cid:80) mk =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; u (cid:17)(cid:33) b φ ( i − m ; j ) ( t − m ; t − t ) . (25) (cid:3) If we set u = 0 in equation (26) we obtain the evolution equation for thetransition probabilities φ ( i − m − ; j ) ( t − m − ; t ): φ ( i − m − ; j ) ( t − m − ; t )= δ i j (cid:32) − H mi (cid:16) m (cid:88) k =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; t − t (cid:17)(cid:33) + (cid:88) i ∈ E t (cid:88) t = t +1 q i i (cid:0) m (cid:88) k =0 i − k − · (cid:0) t − k − t − k − t − t − m − (cid:1) ; t (cid:1) φ ( i − m ; j ) ( t − m ; t − t ) . (26)As it is possible to see from the previous propositions, the transition prob-abilities are function of the last m + 1 states and times of transition, then ourISMC model can also be considered as a type of m + 1 order semi-Markovchain model. Anyway higher order semi-Markov models are associated with adramatic increase in the model parameters that rarely is possible to estimatefrom real data. The ISMC model overcomes this problem because the influ-ence of the past states and times is summarized efficiently with the introduc-tion of the index process that computes a weighted average of the past states.This is particularly evident once we remark that the b φ ( i − m − ; j ) ( t − m − ; u, t )are the unknown functions of formula (16) whereas the remaining quanti-ties are known once we dispose of the indexed semi-Markov kernel Q m ( v ; t )which is the only quantity to be estimated from data. From formula (5) wesee that the influence of past states and times is summarized by the indexprocess V m ( t ). Having estimated Q m ( v ; t ) it is possible, first to computethe sojourn time cumulative distribution H mi ( v ; · ) through formula (11) andthen, by solving equation (16) the transition probability are obtained.
3. Database and commercial wind turbine
As in our previous works [1, 5, 16] we used a free database of wind speedsampled in a weather station situated in Italy at N 45 28’ 14,9” − E 9 22’10ate Wind speed range m/s > Table 1: Wind speed discretization m of altitude. The station uses a combined speed-directionanemometer at 22 m above the ground. It has a measurement range thatgoes from 0 to 60 m/s , a threshold of 0,38 m/s and a resolution of 0,05 m/s .The station processes the speed every 10 minute in a time interval rangingfrom 25/10/2006 to 28/06/2011. During the 10 minutes are performed 31sampling which are then averaged in the time interval. In this work, we usethe sampled data that represents the average of the modulus of the windspeed ( m/s ) without a considered specific direction. This database is thencomposed of about 230thousands wind speed measures ranging from 0 to 16 m/s .To be able to model the wind speed as a semi-Markov process, the statespace of wind speed has been discretized. In the example shown in this workwe discretize wind speed into 8 states (see Table 1) chosen to cover all thewind speed distribution. This choice is done by considering a trade off be-tween accuracy of the description of the wind speed distribution and numberof parameters to be estimated. An increase of the number of states betterdescribes the process but requires a larger dataset to get reliable estimatesand it could also be not necessary for the accuracy needed in forecastingfuture wind speeds. Note also that, in the database used in this work, thereare very few cases where the wind speed exceeds 7 m/s . We stress that thediscretization should be chosen according to the database to be used.We apply our model to a real case of energy production. For this reasonwe choose a commercial wind turbine, a 10 kW Aircon HAWT with a powercurve given in Figure 1. The power curve of a wind turbine represents howmuch energy it produces as a function of the wind speed. In this case, see11 P o w e r ( k W ) wind speed (m/s) Figure 1: Power curve of the 10 kW Aricon HAWT wind turbine.
Figure 1, there is a cut in speed at 2 m/s , instead the wind turbine producesenergy almost linearly from 3 to 10 m/s , then, with increasing wind speedthe energy production remains constant until the cut out speed at 32 m/s ,in which the wind turbine is stopped for structural reason. Then the powercurve acts as a filter for the wind speed. In the database used for our analysisthe wind speed does never exceed 16 m/s and it is seldom over 8 m/s , thisis why the discretization is performed according to Table 1 and the analysisnever reached the cut out speed.Through this power curve we can know how much energy is produced asa function of the wind speed at a given time.
4. Reliability theory for the ISMC model of wind speed
In this section we define and compute reliability measures for the ISMCmodel.Let E be partitioned into sets U and D , so that: E = U ∪ D, ∅ = U ∩ D, U (cid:54) = ∅ , U (cid:54) = E. The subset U contains all ”Up” states in which the system is working andsubset D all ”Down” states in which the system is not working well or hasfailed. In the wind speed model the Up states are those for which the wind12peed is sufficiently high to allow the production of energy or not excessivehigh such that the turbine should be turned off.In the following we present both the typical indicators used in reliabilitytheory and also their application.The three indicators that we evaluate are the availability, reliability andmaintainability functions, they were extensively studied among others by[17, 18] and related to semi-Markov models by [2, 3, 13, 5, 19, 20].(i) the point wise availability function b A with backward time giving theprobability that the system is working on at time t whatever happens on(0 , t ].In our model the availability is defined as follows: b A i − m − ( t − m − ; u, t )= P (cid:2) Z ( t ) ∈ U | T > t + u, ( J , T ) − m − = (( i , t ) − m − ) (cid:3) , (27)and gives the probability that at time t the wind turbine produces energyconditional on the last m + 1 wind speed values i − m − registered at the times t − m − and on the no occurrence of next transition up to time t + u .(ii) the reliability function b R with backward time giving the probabilitythat the system was always working from time t to time t : In our model theavailability is defined as follows: b R i − m − ( t − m − ; u, t )= P (cid:2) Z ( u ) ∈ U, ∀ u ∈ ( t , t ] | T > t + u, ( J , T ) − m − = ( i , t ) − m − (cid:3) (28)and gives the probability that the wind turbine will always produce energyfrom time t up to time t conditional on the last m + 1 wind speed values i − m − registered at the times t − m − and on the no occurrence of next tran-sition up to time t + u .(iii) the maintainability function b M with backward time giving the prob-ability that the system will leave the set D within the time t being in D attime t : b M i − m − ( t − m − ; u, t )= P (cid:2) ∃ u ∈ ( t , t ] : Z ( u ) ∈ U | T > t + u, ( J , T ) − m − = ( i , t ) − m − (cid:3) (29)and gives the probability that the turbine will produce energy at least oncefrom time t up to time t conditional on the last m + 1 wind speed values13 − m − registered at the times t − m − and on the no occurrence of any transi-tion up to time t + u .These three probabilities can be computed in the following way if theprocess is an indexed semi-Markov chain of kernel Q m ( v ; t ).(i) the point wise availability function b A i − m − ( t − m − ; u, t ) :to compute these probabilities it is sufficient to use the following formula: b A i − m − ( t − m − ; u, t ) = (cid:88) j ∈ U b φ ( i − m − ; j ) ( t − m − ; u, t ) . (30)(ii) the reliability function b R m i − m − ( t − m − ; u, t )to compute these probabilities, we will now work with another cumulatedkernel ˆQ m ( v ; t ) for which all the states of the subset D are changed intoabsorbing states by considering the following transformation:ˆ p mi,j ( v ) = i ∈ D, j = i i ∈ D, j / ∈ Dp mi,j ( v ) otherwise. (31) b R i − m − ( t − m − ; u, t ) is given by solving the evolution equation of the indexedsemi-Markov chain but now with the kernel ˆ Q mij ( v ; t ) = ˆ p mij ( v ) G mij ( v ; t ).The related formula will be: b R ( i − m − ; j ) ( t − m − ; u, t ) = (cid:88) j ∈ U b ˆ φ i − m − ( t − m − ; u, t ) (32)where b ˆ φ ( i − m − ; j ) ( t − m − ; u, t ) are the transition probabilities of the processwith all the states in D that are absorbing, i.e. with cumulated kernel ˆ Q .(iii) the maintainability function b M i − m − ( t − m − ; u, t ):to compute these probabilities we will now work with another cumulatedkernel ˜ Q = ( ˜ Q mij ( v ; t )) for which all the states of the subset U are changedinto absorbing states by considering the following transformation:˜ p mij ( v ; t ) = i ∈ U, j = i i ∈ U, j (cid:54) = ip mij ( v ; t ) otherwise. (33)14 M i − m − ( t − m − ; u, t ) is given by solving the evolution equation of an indexedsemi-Markov chain but now with the cumulated kernel ˜ Q mij ( v ; t ) = ˜ p mij ( v ) · G mij ( v ; t ).The related formula for the maintainability function will be: b M i − m − ( t − m − ; u, t ) = (cid:88) j ∈ U b ˜ φ ( i − m − ; j ) ( t − m − ; u, t ) (34)where b ˜ φ ( i − m − ; j ) ( t − m − ; u, t ) are the transition probability of the process withall the states in U that are absorbing, i.e. with cumulated kernel ˜ Q .Unfortunately it is impossible to give a graphical representation of thesethree indicators because they depends on too many parameters ( m + 1 states i − m − , m + 1 times t − m − and one backward value v ), then in order to showthe capacity of the model to correctly reproduce the reliability measures wedefine unconditional reliability measures and we compute them by evaluatingthe discrepancies between real data and model results.Let us define the set H x ,t ( v ) = { ( i − m − , t − m − ) | i = x , t = y , V m ( t ) = v } , then we define the unconditional availability function b A H x ,t ( v ) ( u, t )= P (cid:2) Z ( t ) ∈ U | T > t + u, ( J , T ) − m − = (( a , s ) − m − ) ∈ H x ,t ( v ) (cid:3) . (35)The unconditional availability gives the probability that at time t the windturbine produces energy conditional on the occupancy at the current time t of the state i with the index process having value v and the backward valueequal to u . By using properties of the conditional probabilities it is easy torealize that b A H x ,t ( v ) ( u, t ) = (cid:88) (( a , s ) − m − ) ∈H x ,t ( v ) P [( J , T ) − m − = (( a , s ) − m − )] × P (cid:2) Z ( t ) ∈ U | T > t + u, ( J , T ) − m − = (( a , s ) − m − ) (cid:3) = (cid:88) (( a , s ) − m − ) ∈H x ,t ( v ) µ (( a , s ) − m − ) b A a − m − ( s − m − ; u, t ) . (36)15here µ (( a , s ) − m − ) is the initial distribution of the states a − m − occupiedat times s − m − .The same definitions apply for the reliability and the maintainabilityfunctions.We define the unconditional reliability function b R H x ,t ( v ) ( u, t )= P (cid:2) Z ( u ) ∈ U, ∀ u ∈ ( t , t ] | T > t + u, ( J , T ) − m − = (( a , s ) − m − ) ∈ H x ,t ( v ) (cid:3) . (37)and it results that b R H x ,t ( v ) ( u, t ) = (cid:88) (( a , s ) − m − ) ∈H x ,t ( v ) µ (( a , s ) − m − ) b R a − m − ( s − m − ; u, t ) . (38)The unconditional maintainability function is defined as follows b M H x ,t ( v ) ( u, t )= P (cid:2) ∃ u ∈ ( t , t ] : Z ( u ) ∈ U | T > t + u, ( J , T ) − m − = (( a , s ) − m − ) ∈ H x ,t ( v ) (cid:3) . (39)and it results that b M H x ,t ( v ) ( u, t ) = (cid:88) (( a , s ) − m − ) ∈H x ,t ( v ) µ (( a , s ) − m − ) b M a − m − ( s − m − ; u, t ) . (40)The unconditional reliability measures can be graphically represented be-cause if we set the current time t = 0, then, they depend only on threeparameters ( i , u, v ). In order to verify the validity of our model, we comparethe behaviour of these indicators for real and synthetic data. The indicatorsof the synthetic data are computed averaging over 500 different trajectories.The number of trajectories is chosen to have stable results.The unconditional reliability measures are plotted by varying the backwardtime u , the index process v and maintaining constant the initial state i .The numeric choice of each parameter is given only for graphical reasons, inorder to show the maximum number of curves without overlaps. For all thefigures we use the notation of indicating with a continuous line the indicatorscomputed from real data and with a dashed line those computed from simu-lated data. As numeric indicator to compare the gap between the curves wecompute the root mean square error (RMSE) between the indicator applied16
0 60 1000.4 0.6 0.8 Time (Minutes) P r obab ili t y i=3 v=3 u=1i=3 v=3 u=3i=3 v=2 u=1i=3 v=2 u=3 Figure 2: Comparison of unconditional availability functions for real and simulated data to the real data and the 500 simulated trajectories. The mean square erroris defined as follows: b RM SE i,v = (cid:118)(cid:117)(cid:117)(cid:116) (cid:88) h =1 (cid:0) b I reali,v − b I hi,v (cid:1) (41)where b I reali,v stands for the indicators estimated form real data and b I hi,v theindicators estimated from each synthetic trajectory.In Figure 2 the unconditional availability functions of real and synthetic10 min min min RM SE , RM SE , RM SE , RM SE , Table 2: Mean square error between the curves of the availability applied to real andsynthetic data
0 60 1000 0.4 0.8 Time (Minutes) P r obab ili t y i=3 v=2 u=1i=3 v=2 u=2i=3 v=4 u=1i=3 v=4 u=3 Figure 3: Comparison of unconditional reliability functions for real and simulated data data are compared. The comparison is made by varying the backward time u and the index process v and maintaining constant the starting state i . InTable 2 there are the RMSE values computed between real and simulatedcurves of Figure 2. As it is possible to note there is a slight increasing trendof the RMSE at the increasing of time. This consideration can be extendedalso to the unconditional reliability and maintainability functions.Figure 3 shows the reliability function for real data compared with thosesimulated. The plotting procedure is the same as for the previous figure.The theoretical trend of the reliability function is confirmed, the probability10 min min min RM SE , RM SE , RM SE , RM SE , Table 3: Mean square error between the curves of the reliability applied to real andsynthetic data
0 60 1000.2 0.6 1 Time (Minutes) P r obab ili t y i=2 v=2 u=1i=2 v=2 u=3i=2 v=4 u=1i=2 v=4 u=3 Figure 4: Comparison of unconditional maintainability functions for real and simulateddata decreases at the increasing of the time interval. In Table 3, instead there arethe RMSE values computed between real and simulated curves of Figure 3.The maintainability function is plotted in Figure 4. As the previous fig-ures, this one shows the comparison of the maintainability for real and simu-lated data varying the backward time and the index process and maintainingconstant the starting state. Table 4 shows the RMSE values computed be-tween real and simulated curves of Figure 4.It is possible to note that all the indicators plotted above (availability,10 min min min RM SE , RM SE , RM SE , RM SE , Table 4: Mean square error between the curves of the maintainability applied to real andsynthetic data
Table 5: Root mean square error percentage between real and synthetic data for availabil-ity, reliability and maintainability function computed with different models reliability and maintainability) depend strongly on the backward time andon the index process. In fact, all the probabilities have different values alsoif only the index process v is changed keeping constant the backward time b and the starting state i . For example, from Figure 3 it is possible to see thatfor all s b R H , (4) (1 , s ) > b R H , (2) (1 , s ) ∀ s ∈ [0 , , and in particular b R H , (4) (1 ,
40) = 0 ,
77 and b R H , (2) (1 ,
40) = 0 ,
41. Thisreveals that it is important to dispose of a model that is able to distinguishbetween these different situations which are determined only from a differ-ent duration of permanence in the initial state i and different values of theaverage of past m + 1 visited states. Models based on Markov chains or clas-sical semi-Markov chain are unable to capture this important effect that ourindexed semi-Markov chain reproduces accordingly to the real data. To mo-tivate well this statement, in Table 5 we show the percentage RMSE betweenall the indicators (availability, reliability and maintainability) evaluated forreal and synthetic data computed through different models. Particularly weuse, for this comparison, a first order Markov chain, a first order semi-Markovchain, a second order in state and duration semi-Markov chain and the ISMCmodel. From the Table it is clear that the best model is the ISMC model.
5. Conclusion
In our previous work, we presented a new stochastic model for the gen-eration of synthetic wind speed data based on a semi-Markov approach butincluding a new and important random variable that able us to capture wellthe behaviour of the wind speed. Here we compute, for the first time, typi-cal indicators in reliability theory for wind speed phenomenon by using theISMC model. 20n order to check the validity of the presented model, we have comparedthe behaviour of these indicators for real and simulated data. To do this,we applied our model to a real case of energy production filtering real andsimulated data by the power curve of a commercial wind turbine. The resultsshow that the proposed model is able to reproduce the behaviour of real databy exhibiting the dependence of the reliability indicators on the backwardtime and the index process.The indications provided by the model are of importance for assessingthe suitability of a location for the wind farm installation as well as for theplanning of a preventive maintenance policy. We have also shown that theISMC model reproduces the behavior of real data much better than a firstorder Markov chain, a first order semi-Markov chain and a second order instate and duration semi-Markov chain.
References [1] G. D’Amico, F. Petroni, and F. Prattico, Wind speed modeled as anindexed semi-Markov process,
Environmetrics (2013) 367-376.[2] V. Barbu, M. Boussemart and N. Limnios, Discrete Time Semi-MarkovModel for Reliability and Survival Analysis, Communications in Statis-tics, Theory and Methods (2004) 2833-2868.[3] A. Blasi, J. Janssen, R. Manca, Numerical treatment of homogeneousand non-homogeneous semi-Markov reliability models. Communicationsin Statistics, Theory and Methods (2004) 697-714.[4] N. Limnios, G. Opri¸san, An introduction to Semi-Markov Processes withApplication to Reliability, In D.N. Shanbhag and C.R. Rao, eds., Hand-book of Statistics (2003) 515-556.[5] G. D’Amico, F. Petroni, and F. Prattico, Reliability Measures of Second-Order Semi-Markov Chain Applied to Wind Energy Production, Journalof Renewable Energy (2013) Article ID 368940.[6] G. D’Amico, J. Janssen, and R. Manca, Semi-Markov Reliability Mod-els with Recurrence Times and Credit Rating Applications,
Journal ofApplied Mathematics and Decision Sciences
Article ID 625712, (2009)17 pages. 217] G. D’Amico, F. Petroni, A semi-Markov model with memory for pricechanges,
Journal of Statistical Mechanics
P12009 (2011).[8] G. D’Amico, Age-usage semi-Markov models,
Applied MathematicaModelling , (2011) 4354-4366.[9] J. Janssen and R. Manca, Applied semi-Markov processes, Springer ,New York, 2006.[10] V. Barbu and N. Limnios. Semi-Markov Chains and Hidden Semi-Markov Models Toward Applications. Springer-Verlag New York Inc,2008.[11] G. D’Amico, F. Petroni, A semi-Markov model for price returns.
PhysicaA (2012) 4867-4876.[12] G. Ciardo, R. Marie, B. Sericola, K.S. Trivedi, Performability Analysisusing Semi-Markov Reward Processes,
IEEE Transactions on Comput-ers , C − (10) (1990) 1251-1264.[13] G. D’Amico, The crossing barrier of a non-homogeneous semi-Markovchain, Stochastics: An International Journal of Probability and Stochas-tic Processes , (6) (2009) 589-600.[14] A. Csenki, Transition analysis of semi-Markov reliability models - a tu-torial review with emphasis on discrete-parameter approaches, In S. Os-aki, editor, Stochastic Models in Reliability and Maintenance, (2002)219-251, Springer, Berlin.[15] G. D’Amico, J. Janssen, and R. Manca, Duration Dependent Semi-Markov Models,
Applied Mathematical Sciences (42) (2011) 2097-2108.[16] G. D’Amico, F. Petroni and F. Prattico, First and second order semi-Markov chains for wind speed modeling, Physica A (2013) 1194-1201.[17] R.E. Barlow, F. Prochan,
Statistical Theory of Reliability and Life Test-ing: Probabilistic Models.
Holt, Rinehart and Winston, New York, 1975.[18] G. Lisnianski, G. Levitin,
Multi-State System Reliability.
World Scien-tific, N.J. (2003). 2219] R. Perez-Ocon, I. Torres-Castro, A reliability semi-Markov Model in-volving geometric processes,
Applied Stochastic Models in Business andIndustry (2) (2002) 157-170.[20] N. Limnios, G. Opri¸san, A unified approach for reliability and performa-bility evaluation of semi-Markov systems, Applied Stochastic Models inBusiness and Industry15