Transcript Chapter 10 Correlation Analysis SIMULATION MODELING AND ANALYSIS WITH ARENA

SIMULATION MODELING AND ANALYSIS
WITH ARENA
T. Altiok and B. Melamed
Chapter 10
Correlation Analysis
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
1
What is Correlation Analysis?
• Correlation Analysis is a modeling and analysis
approach that straddles both Input Analysis and
Output Analysis
• Correlation Analysis consists of two activities:
1. modeling of correlated stochastic processes
2. studying the impact of correlations on performance
measures of interest via Sensitivity Analysis
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
2
Correlation in Input Analysis
• Correlation Analysis as part of Input Analysis is simply an
approach to modeling and data fitting that
• insists on high-quality models incorporating temporal dependence
• strives to fit correlation-related statistics in a systematic way
• To set the scene, consider a stationary time series , {X n }n¥ = 0,
that is, all statistics remain unchanged under the passage of time
• in particular, all X n share a common mean, m X and common variance, s X2
• to fix the ideas, suppose that {X n } is to be used to model inter-arrival times
at a queue (in which case the time series is non-negative)
• What statistical aspects of
{X n } n¥ = 0
should be carefully modeled?
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
3
Statistical Signatures
• Let a collection of statistics of a random variable or time series
be referred to as a statistical signature (signature, for short)
• it is often possible to order signatures by “strength”, for example,
signatures obviously become stronger under inclusion
• To clarify the signature strength notion, consider a time series of
inter-arrival times, X = {X n } , and the following set of signatures
in increasing strength
• the mean, m X , is a “minimal” signature, since its reciprocal, l = 1/ m ,
X
is the arrival rate – a key statistic in queueing models
• the mean, m X , and variance, s X2 , is a stronger signature
• adding moments of the inter-arrival distribution, such as the skewness
and kurtosis, yields an even stronger signature
• the (marginal) distribution, FX , determines all its moments,
and so is stronger than all of the above
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
4
A Very Strong Signature
• Given a stationary empirical time series, our goal here is to to fit
a particular very strong signature to a time series, X = {X n ,}
which includes both of the following statistics:
• the marginal distribution, F
X
E éêX n X n + t
• the autocorrelation function, r X (t ) = ë
2
• The marginal distribution,
s
ù- m 2
ú
X
û
, t = 1, 2 , …
X
F
X
• is a first-order statistic of { X n }, that is, it involves only a single
random variable from { X n } (by stationarity)
• is estimated by an empirical histogram, Hˆ
• The autocorrelation function, r X (t )
Y
• is a second-order statistic of { X n } , that is, it involves pairs of lagged
random variables from { X n } , and
• serves as a statistical proxy for temporal dependence in { X n } , where each
correlation coefficient, r X (t ) = r [X n , X n + t ] , measures of linear dependence
• is estimated by some estimator rˆ X (t ), t = 1, ...,T
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
5
Correlation in Output Analysis
• Correlation Analysis as part of Output Analysis is the study of the
sensitivity of output statistics to correlations in model components
• autocorrelation can have a major impact on performance measures
• consequently, they cannot always be ignored merely for the sake of
simplified models
• however, modelers routinely ignore correlations to simplify model
construction and its analysis
• A motivating example from the domain of queueing systems
will illustrate the peril of ignoring correlations uncritically
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
6
Example: correlation Impact
• Consider a workstation operating as an M/M/1 system with job
arrival rate l and processing (service) rate m, such that m> l
• since all job arrivals and processing times are mutually independent,
all corresponding autocorrelations and cross-correlations are identically zero
• the system is stable with utilization u = l / m < 1
• it is known that the equilibrium mean flow time is E[S
] = 1/ (m- l )
M /M /1
and so the mean waiting time in the buffer is E[W M / M / 1 ] = 1/ (m- l ) - 1/m
• Next, modify the arrival process from a Poisson process to a
(possibly autocorrelated) TES process, yielding a TES/M/1 system
• The merit of TES processes is that they simultaneously admit arbitrary
marginal distributions and a variety of autocorrelation functions
• in particular, we can select TES inter-arrival processes with the same interarrival time distribution as in the Poisson process (i.e., exponential with rate
parameter l ), but with autocorrelated inter-arrival times, yielding some
TES/M/1 equilibrium mean waiting time in the buffer, E[W
]
T ES /M /1
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
7
Example: correlation Impact (Cont.)
• We wish to gauge the impact of autocorrelations in the job arrival
stream on mean waiting times via the relative deviation
E[W
d (r (1)) =
T ES /M /1
E[W
] - E[W
M /M /1
M /M /1
]
]
• the relative deviation is viewed as a function of the lag-1 autocorrelation,
r (1) , in the TES arrival process
• The table below displays the relative deviations for two
representative cases:
• l = 0.25 and m = 1 (light traffic regime with utilization u = 0.25)
• l = 0.80 and m = 1 (heavy-traffic regime with utilization u = 0.80 )
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
8
Introduction to TES Modeling
• The definition of a TES process involves two related
stochastic processes:
• an auxiliary process, called the background process
• a target process, called the foreground process
• The two processes operate in lockstep in the sense that
they are connected by a deterministic transformation
• more specifically, the state of the background process is mapped to a state
of the foreground process
• this is done in such a way that the foreground process has a prescribed
marginal distribution and a prescribed autocorrelation function
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
9
Modulo-1 Arithmetic
• The definition of TES processes makes use of a simple
mathematical operation, called modulo-1 arithmetic
• modulo-1 arithmetic is arithmetic restricted to the familiar fractions
(values in the interval [0,1), with the value 1 excluded)
• the notation x = x - max{integer n : n £ x } is used to denote the
fractional value of any number
• note that fractional values are defined for any real number
(positive as well as negative)
• Examples:
• for zero, we simply have 0  0
• for positive numbers, we have the familiar fractional values, for example,
0.6 = 1.6 = 2.6 = ... = 0.6
• for negative values, the fractional part is the complementary value relative
to one, for example,
.
0.6  1.6  2.6  ... 1 0.6  0.4
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
10
Outline of TES Processes Theory
• The Lemma of Iterated Uniformity is the foundation of
the theory of TES processes
• let U be a uniform random variable on [0,1) and let V be any random
variable independent of U
• then U + V is also uniform on [0,1), regardless of the distribution of V !
• Define a stochastic process
U 0, áU 0 + V 1ñ, áU 0 + V 1 + V 2 ñ, ..., áU 0 + V 1 + ... + V n ñ, ...
• by the Lemma of Iterated Uniformity, each random variable above is
uniform on [0,1)
• furthermore, each could be further transformed into a foreground process
F
- 1
(U 0 ), F
- 1
(áU 0 + V 1 ñ), F
- 1
(áU 0 + V 1 + V 2 ñ), ..., F - 1(áU 0 + V 1 + ... + V n ñ), ...
and by the Inverse Transform Method, each random variable above
will have the prescribed distribution F !
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
11
Background TES Processes
• Define the following random variables:
• let U 0 be a random variable with a uniform distribution on [0,1)
• let {V n } n¥ = 1 be an innovation sequence (that is, any iid sequence of
random variables, independent of U 0 )
• TES background processes come in two flavors:
• a background TES+ process, {U n+ }, is defined by the recursive scheme
ìïU 0,
U n = ïí
ïï áU n+ î
n = 0
+
1 + V n ñ,
n > 0
• a background TES- process, {U n- }, is defined by
ìïU + ,
ï n
Un = í
ïï 1- U + ,
n
ïî
n even
n odd
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
12
Visualizing Background TES Processes
• Background TES processes can be visualized as a random walk
on the unit circle
• Consider a basic TES process, where the innovation variate is uniform on
an interval [L, R ), so its density is a single step of length not exceeding 1
innovation
step
density
innovation
density
origin
+
U
áU n + L ñ
n+1
+
U n+
áU n+ + R ñ
unit circle
origin
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
13
Basic TES Processes
• The following list summarizes qualitatively the effect of the
parameters Land R on the autocorrelation of a basic background
TES process:
• the width, R - L , of the innovation-density support (the region over which
the density is positive) has a major effect on the magnitude of the
autocorrelations: the larger the support, the smaller the magnitude
(in fact, when R - L = 1 , then the autocorrelations vanish altogether)
• the location of the innovation-density support affects the shape of the
autocorrelation function: when the support is not symmetric about the
origin, then the autocorrelation function assumes an oscillating form,
and otherwise it is monotone decreasing
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
14
Basic TES Processes (Cont.)
Autocorrelation function of a basic TES+ process (L, R ) = (- 0.05, 0.05)
(symmetric innovation density and narrow support)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
15
Basic TES Processes (Cont.)
Autocorrelation function of a basic TES+ process (L, R ) = (- 0.2, 0.2)
(symmetric innovation density and wider support)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
16
Basic TES Processes (Cont.)
Autocorrelation function of a basic TES+ process (L, R ) = (- 0.01, 0.05)
(non-symmetric innovation density)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
17
Basic TES Processes (Cont.)
Autocorrelation function of a basic TES- process (L, R ) = (- 0.05, 0.05)
(symmetric innovation density)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
18
Basic TES Processes (Cont.)
Autocorrelation function of a basic TES- process (L, R ) = (- 0.01, 0.05)
(non-symmetric innovation density)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
19
Stitching Transformations
• A background TES process can produce marked visual
“discontinuities” in its sample paths, which are noticeable when
• the innovation density has a narrow support
• successive background variates on the unit circle straddle the circle’s origin
• in the figure below we have a “sudden” drop from a relatively high value to
a relatively small one as the process crosses the origin counter clock-wise
Sample path of a basic TES+ background process with (L, R ) = (- 0.04, 0.055)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
20
Stitching Transformations (Cont.)
• For modeling purposes, we would like sometimes to “smooth”
(“stitch together”) such marked visual “discontinuities”
• To this end, define a family of so-called stitching transformations
ìï u ,
if 0 £ u £ x
ïï
x
ï
S (u ) = í
x
ïï 1 - u
ïï 1 - x , if x £ u £ 1
î
where x is a so-called stitching parameter in the interval [0,1]
• a stitching transformation preserves uniformity, that is,
if U ~ Unif( 0, 1), then S (U ) ~ Unif( 0, 1) for any x Î [0, 1]
x
• therefore, any stitched background TES sequence is also a TES background
sequence (and thus uniformly distributed), albeit a “smoothed” one
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
21
Stitching Transformations (Cont.)
• The graph below displays typical stitching transformations for the
following stitching parameters, x
• for x = 0 , S 0(u ) = 1 - u
• for 0 < x < 1 , S (u ) has a triangular shape
x
• for x = 1 , S 1(u ) = u is the identity
S (u )
1
x
S 1(u ) = u
S 0(u ) = 1 - u
0
1
x
Stitching transformations for S 0 (dashed curve), S 1 (dotted curve) and S
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
x
(solid curve)
22
Stitching Transformations (Cont.)
• The graphs below illustrate the smoothing effect of stitching
Sample path of a basic TES+ background process with (L, R ) = (- 0.04, 0.055)
and without stitching ( x = 1)
Sample path of a basic TES+ background process with (L, R ) = (- 0.04, 0.055)
and with stitching (x = 0.5 )
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
23
Foreground TES Processes
• A foreground TES process is obtained from a background TES
process by a deterministic transformation, D , called a distortion
• a foreground TES+ process, { X n+ }, is of the form X n+ = D (U n+ )
• a foreground TES- process, { X n- }, is of the form X n- = D (U n- )
• In practice, one often applies a stitching transformation followed
by an application of the Inverse Transform method via a
distortion of the form
D (u ) = F
where
-
1
(S x (u ))
• S x is a stitching transformation (often x = 0.5 )
• F is a cdf (typically, F = Hˆ h is an empirical histogram of data vector, h )
• for example, for the exponential cdf, Fexp (x ) = 1- e - l x , the inverse is
- 1
1
Fexp (y ) =
ln(1 - y ), and the stitching transformation might be S
l
0.5
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
24
Foreground TES Processes (Cont.)
• Example: applying the exponential Inverse Transform formula
above to basic background TES+ processes to obtain
foreground TES+ processes
• two basic TES+ background processes are used with, respectively,
(L, R ) = (- 0.05, 0.05) and (L, R ) = (- 0.5, 0.5)
• the Inverse Transform applied to these TES+ background processes
uses the same parameter, l = 1
• The results are shown in the next few foils:
• both foreground TES+ processes have the same exponential marginal
distribution of rate 1, as attested by their histograms
• in contrast, the first foreground process exhibits significant
autocorrelations, while the second has zero autocorrelations,
as a consequence of its iid property (see the corresponding correlograms)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
25
Foreground TES Processes (Cont.)
Sample path (top), histogram (middle) and correlogram (bottom)
of an exponential basic TES+ process with background parameters (L, R ) = (- 0.05, 0.05)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
26
Foreground TES Processes (Cont.)
Sample path (top), histogram (middle) and correlogram (bottom)
of an exponential basic TES+ process with background parameters (L, R ) = (- 0.5, 0.5)
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
27
Generation of TES Sequences
• TES processes are readily generated on a computer via
algorithms that utilize random number generators
(RNG)
• we assume that the availability of a function, called mod1(x),
which implements modulo-1 reduction of any real number,
and returns the corresponding fraction
• For convenience, we separate the generation of TES+
processes from that of TES- processes
• the corresponding algorithms have considerable overlaps
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
28
Generation of
+
TES
Sequences
• Inputs:
• an innovation density fV
• a stitching parameter x
• an inverse distribution F -
• Outputs:
1
// modeler choice
// modeler choice
// often inverse histogram (step) cdf, Ĥ h- 1
• a background TES+ sequence, {U n+ } ¥
n= 0
• a foreground TES+ sequence, { X n+ } ¥
• Algorithm:
n= 0
1. sample a value U 0 , uniform on [0,1), // initial background variate
and set U 0+ ¬ U 0 and n ¬ 0
// more initializations
2. go to Step 6.
// go to generate initial foreground variate
3. set n ¬ n + 1
// bump up running index for next iteration
4. sample a value V from f // sample an innovation variate
V
+
+
5. set U n ¬ mod1(U + V ) // compute a TES+ background variate
n- 1
+
6. set S ¬ S (U n )
// apply a stitching transformation
x
7. set X n+ ¬ F - 1(S )
// compute a TES+ foreground variate
8. go to Step 3.
// loop to generate the next TES+ variate
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
29
Generation of TES Sequences
• Inputs:
• an innovation density f
V
• a stitching parameter x
• an inverse distribution F -
• Outputs:
1
// modeler choice
// modeler choice
// often inverse histogram (step) cdf, Ĥ h- 1
• a background TES- sequence, {U n- } n¥ = 0
• a foreground TES- sequence, { X n- } ¥
• Algorithm:
n= 0
1. sample a value U 0 , uniform on [0,1), // initial background variate
and set U 0- ¬ U 0+ ¬ U 0 and n ¬ 0
// more initializations
2. go to Step 7.
// go to generate initial foreground variate
3. set n ¬ n + 1
// bump up running index for next iteration
4. sample a valueV from fV // sample an innovation variate
5. set U n+ ¬ mod1(U n+ - 1 + V ) // compute a TES+ background variate
6. if n is even, then set U n- ¬ U n+ , // compute a TES- background variate
if n is odd, then set U n- ¬ 1 - U n+
7. set S ¬ S x (U n- )
// apply a stitching transformation
8. set X n- ¬ F - 1(S )
// compute a TES- foreground variate
9. go to Step 3.
// loop to generate the next TES- variate
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
30
Generation of TES Sequences in Arena
• Arena implementation of the algorithm to generate basic TES+
sequences with an exponential distribution (TES- is similar)
assumes that the following parameters are given as inputs:
• a pair of parameters some L and R such that 0 £ L < R < 1 ,
which determine a basic innovation density
ìï
1
,
ïï
ï
R
L
f (x ) = í
V
ïï
ïïî 0,
L£ x< R
otherwise
• a stitching parameter 0 £ x £ 1 (0.5 is typical)
• an inverse of an exponential distribution function,
F - 1 (u ) = -
1
ln(1 - u )
l
for some l > 0
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
31
Arena Model for Basic TES+
Arena model implementing the generation of
basic TES+ sequences with exponential marginal distribution
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
32
Arena Model for Basic
+
TES (Cont.)
Arena Variable module for implementing the generation of
basic TES+ sequence with an exponential marginal distribution
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
33
Arena Model for Basic
+
TES (Cont.)
• Arena variables in the model:
• the variables L and R hold the parameters of the basic innovation density
• the variable xi holds the stitching parameter
• the variable lambda holds the rate parameter of the requisite exponential
distribution
• the variable N holds the running index in the TES+ sequence (initially 0)
• the variable V_N holds an innovation variate
• the variable U_N holds an unstitched TES+ background variate
• the variable US_N holds a stitched TES+ variate
• the variable UP_N holds a stitched TES+ background variate
• the variable X_N holds a TES+ foreground variate
• The Arena Variable module
• lists all model parameters and variables and their initial values, if any
• those requiring initialization are identified by a 1 rows button label
under the Initial Values heading), for example UP_N is initialized to a
value between 0 and 1
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
34
Correlation Analysis Example
• Consider a workstation subject to mutually independent failures
• in the case, the workstation can be modeled as an M/G/1 queueing system,
where the processing time is, in fact, the process completion time,
consisting of all the failures experienced by a job on the machine
• The mean job waiting time is given by the modified P-K formula
2
E[W M / G / 1 ] =
where
l c (g C2 + 1)
2 (1 - l c )
• l is the arrival rate
• c = x (1 + d y ) is the average process completion time
•
g C2
2
2
= x d y / c is the squared coefficient of variation of the process
2
completion time, where y is the second moment of repair times
• The probability that the machine is occupied (processing or down)
is given by, Pr{B } = l c and for stability we assume Pr{B } < 1
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
35
Correlation Analysis Example (Cont.)
• The table below displays the relative deviations
d ( r (1)) =
as function of
Pr{ B }
E[W T ES / M / 1 ] - E[W M / G / 1 ]
E[W M / G / 1 ]
, where
• Pr{B } = 0.66 for x = 0.25
x = 0.31
• Pr{B } = 0.81 for
Pr{B }
0.66
0.81
0.00
15.8
36.3
Lag-1 Autocorrelation of Time-to-Failure
r (1)
0.14
0.36
0.48
0.60
0.68
0.74
0.99
11.4% 51.3% 90.1% 260%
543.7% 3,232% 4,115%
11%
151.3% 229.8% 564.7% 766%
3,429% 7,800%
Relative deviations of mean waiting time in a workstation with failures/repairs
Altiok / Melamed Simulation Modeling and Analysis with Arena
Chapter 10
36