1
2# This library is free software; you can redistribute it and/or
3# modify it under the terms of the GNU Library General Public
4# License as published by the Free Software Foundation; either
5# version 2 of the License, or (at your option) any later version.
6#
7# This library is distributed in the hope that it will be useful,
8# but WITHOUT ANY WARRANTY; without even the implied warranty of
9# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10# GNU Library General Public License for more details.
11#
12# You should have received A copy of the GNU Library General
13# Public License along with this library; if not, write to the
14# Free Foundation, Inc., 59 Temple Place, Suite 330, Boston,
15# MA  02111-1307  USA
16
17
18###############################################################################
19# FUNCTION:             STYLIZED FACTS GUI:
20#  .stylizedFactsGUI     Opens a GUI for stylized facts
21#  .lmacfPlot            Estimates and displays the long memory ACF
22#  .logpdfPlot           Displays a pdf plot on logarithmic scale(s)
23#  .ccfPlot              Displays tailored cross correlation function plot
24#  .qqgaussPlot          Displays a tailored Gaussian quantile-quantile plot
25###############################################################################
26
27
28.lmacfPlot <-
29function(x, lag.max = max(2, floor(10*log10(length(x)))),
30    ci = 0.95, type = c("both", "acf", "hurst"), labels = TRUE,
31    trace = TRUE, ...)
32{
33    # A function implemented by Diethelm Wuertz
34
35    # Description:
36    #   Evaluate and display long memory autocorrelation Function.
37
38    # Arguments:
39    #   x - an uni- or multivariate return series of class 'timeSeries'
40    #       or any other object which can be transformed by the function
41    #       'as.timeSeries()' into an object of class 'timeSeries'.
42    #   labels - a logical flag, by default true. Should a default
43    #       main title and labels addet to the plot?
44
45    # FUNCTION:
46
47    # Settings:
48    if (!is.timeSeries(x)) x = as.timeSeries(x)
49    Units = colnames(x)
50    X = x
51
52    # Match Arguments:
53    type = match.arg(type)
54
55    # Labels:
56    if (labels) {
57        main1 = ""
58        xlab1 = "lag"
59        ylab1 = "ACF"
60        main2 = ""
61        xlab2 = "log lag"
62        ylab2 = "log ACF"
63    } else {
64        main1 = xlab1 = ylab1 = ""
65        main2 = xlab2 = ylab2 = ""
66    }
67
68    Fitted = list()
69    Hurst = NULL
70    DIM = dim(X)[2]
71    for (i in 1:DIM) {
72
73        # Get Data:
74        x.ret = as.matrix(X)[, i]
75        x = abs(x.ret)
76        if (labels) main1 = main2 = Units[i]
77
78        # ACF:
79        z = acf(x, lag.max = lag.max, type = "correlation", plot = FALSE)
80        z$acf[1] = 0
81        cl = qnorm(0.5 + ci/2)/sqrt(z$n.used)
82        z.min = min(z$acf, -cl)
83
84        # lin-lin plot excluding one:
85        x = seq(0, lag.max, by = 1)
86        y = z$acf
87        if (type == "both" | type == "acf") {
88            plot(x = x[-1], y = y[-1], type = "l", main = main1,
89                col = "steelblue", xlab = xlab1, ylab = ylab1,
90                xlim = c(0, lag.max), ylim = c(-2*cl, max(y[-1])), ...)
91            # abline(h = 0, lty = 3)
92            if (trace) {
93                cat ("\nLong Memory Autocorrelation Function:")
94                    paste (cat ("\n  Maximum Lag        "), cat(lag.max))
95                    paste (cat ("\n  Cut-Off ConfLevel  "), cat(cl))
96            }
97            ACF = acf(x.ret, lag.max = lag.max, plot = FALSE)$acf[,,1]
98            lines(x = 1:lag.max, y = ACF[-1], type = "l", col = "steelblue")
99            lines(x = c(-0.1, 1.1)*lag.max, y = c(+cl, +cl), lty = 3,
100                col = "darkgrey")
101            lines(x = c(-0.1, 1.1)*lag.max, y = c(-cl, -cl), lty = 3,
102                col = "darkgrey")
103        }
104
105        # log-log Plot of ACF:
106        x = x[y > cl]
107        y = y[y > cl]
108        # log-log:
109        if (length(x) < 10) {
110            Fit = c(NA, NA)
111            hurst = NA
112            cat("\n  The time series exhibits no long memory! \n")
113        } else {
114            Fit = lsfit(log(x), log(y))
115            fit = unlist(Fit)[1:2]
116            hurst = 1 + fit[2]/2
117            if (type == "both" | type == "hurst") {
118                plot(x = log(x), y = log(y), type = "l", xlab = xlab2,
119                    ylab = ylab2, main = main2, col = "steelblue", ...)
120                Text = paste("Hurst Exponent:", signif(hurst, 3))
121                mtext(Text, side = 4, adj = 0, col = "darkgrey", cex = 0.7)
122                # Grid:
123                if (labels) grid()
124            }
125            ### fit = l1fit(log(x), log(y))$coefficients
126            abline(fit[1], fit[2], col = 1)
127            if (trace) {
128                paste (cat ('\n  Plot-Intercept     '), cat(fit[1]))
129                paste (cat ('\n  Plot-Slope         '), cat(fit[2]))
130                paste (cat ('\n  Hurst Exponent     '), cat(hurst), cat("\n"))
131            }
132        }
133
134        # Save Results:
135        if (DIM == 1) Fitted = Fit else Fitted[[i]] = Fit
136        Hurst = c(Hurst, hurst)
137    }
138
139    # Return Value:
140    invisible(list(fit = Fitted, hurst = Hurst))
141}
142
143
144# ------------------------------------------------------------------------------
145
146
147.logpdfPlot =
148function(x, breaks = "FD", type = c("lin-log", "log-log"),
149doplot = TRUE, labels = TRUE, ...)
150{   # A function implemented by Diethelm Wuertz
151
152    # Description:
153    #   Displays a pdf plot on logarithmic scale(s)
154
155    # Details:
156    #   Returns a pdf plot on a lin-log scale in
157    #   comparison to a Gaussian density plot
158    #   and return the break-midpoints and the
159    #   counts obtained from the histogram of
160    #   the empirical data.
161
162    # Arguments:
163    #   x - an uni- or multivariate return series of class 'timeSeries'
164    #       or any other object which can be transformed by the function
165    #       'as.timeSeries()' into an object of class 'timeSeries'.
166    #   labels - a logical flag, by default true. Should a default
167    #       main title and labels addet to the plot?
168
169    # FUNCTION:
170
171    # Settings:
172    if (!is.timeSeries(x)) x = as.timeSeries(x)
173    Units = colnames(x)
174
175    # Select Type:
176    type = match.arg(type)
177
178    # Labels:
179    if (labels) {
180        if (type == "lin-log") {
181            main = "log PDF"
182            xlab = "x"
183            ylab = "log PDF"
184        } else if (type == "log-log") {
185            main = "log PDF"
186            xlab = "log x"
187            ylab = "log PDF"
188        }
189    } else {
190        main = xlab = ylab = ""
191    }
192
193    X = x
194    DIM = ncol(X)
195
196    for (i in 1:DIM) {
197
198        # Get Data:
199        x = as.vector(X[, i])
200        if (labels) main = Units[i]
201
202        # Lin-Log Plot:
203        if (type == "lin-log") {
204
205            # Histogram PDF:
206            result = hist(x, breaks = breaks, plot = FALSE)
207            prob.counts = result$counts/sum(result$counts) /
208                diff(result$breaks)[1]
209            histogram = list(breaks = result$breaks, counts = prob.counts)
210
211            # Histogram Count & Break-Midpoints:
212            yh = histogram$counts
213            xh = histogram$breaks
214            xh = xh[1:(length(xh)-1)] + diff(xh)/2
215            xh = xh[yh > 0]
216            yh = log(yh[yh > 0])
217            if (doplot) {
218                par(err = -1)
219                plot(xh, yh, type = "p", pch = 19, col = "steelblue",
220                    main = main, xlab = xlab, ylab = ylab, ...)
221                Text = "Scales: log-log"
222                mtext(Text, side = 4, adj =0, col = "darkgrey", cex = 0.7)
223            }
224
225            # Compare with a Gaussian Plot:
226            xg = seq(from = xh[1], to = xh[length(xh)], length = 301)
227            yg = log(dnorm(xg, mean(x), sqrt(var(x))))
228            if (doplot) {
229                par(err = -1)
230                lines(xg, yg, col = "brown")
231            }
232
233            # Return Value:
234            result = list(breaks = xh, counts = yh, fbreaks = xg,
235                fcounts = yg)
236        }
237
238        # Log-Log Plot:
239        if (type == "log-log") {
240
241            # Histogram PDF:
242            result = hist(x, breaks = breaks, plot = FALSE)
243            prob.counts = result$counts/sum(result$counts) / diff(result$breaks)[1]
244            histogram = list(breaks = result$breaks, counts = prob.counts)
245
246            # Histogram Count & Breaks:
247            yh = histogram$counts
248            xh = histogram$breaks
249            xh = xh[1:(length(xh)-1)] + diff(xh)/2
250            xh = xh[yh > 0]
251            yh = yh[yh > 0]
252            yh1 = yh[xh < 0]
253            xh1 = abs(xh[xh < 0])
254            yh2 = yh[xh > 0]
255            xh2 = xh[xh > 0]
256            if (doplot) {
257                plot(log(xh1), log(yh1), type = "p", pch = 19,
258                    col = "darkgreen",
259                    main = main, xlab = xlab, ylab = ylab, ...)
260                Text = "Scales: log-log"
261                mtext(Text, side = 4, adj =0, col = "darkgrey", cex = 0.7)
262                par(err = -1)
263                points(log(xh2), log(yh2), col = 2, ...)
264            }
265
266            # Compare with a Gaussian plot:
267            xg = seq(from = min(xh1[1], xh[2]),
268                to = max(xh1[length(xh1)], xh2[length(xh2)]), length = 301)
269            xg = xg[xg > 0]
270            yg = log(dnorm(xg, mean(x), sqrt(var(x))))
271            if (doplot) {
272                par(err = -1)
273                lines(log(xg), yg, col = "brown")
274            }
275
276            # Result:
277            result = list(breaks = c(xh1, xh2), counts = c(yh1, yh2),
278                fbreaks = c(-rev(xg), xg), fcounts = c(-rev(yg), yg))
279        }
280
281        # Grid:
282        if (labels) grid()
283
284    }
285
286    # Return Value:
287    invisible(result)
288}
289
290
291# ------------------------------------------------------------------------------
292
293
294.qqgaussPlot <-
295function(x, span = 5, col = "steelblue", labels = TRUE, ...)
296{
297    # A function implemented by Diethelm Wuertz
298
299    # Description:
300    #   Returns a simple Quantile-Quantile plot.
301
302    # Arguments:
303    #   x - an uni- or multivariate return series of class 'timeSeries'
304    #       or any other object which can be transformed by the function
305    #       'as.timeSeries()' into an object of class 'timeSeries'.
306    #   labels - a logical flag, by default true. Should a default
307    #       main title and labels addet to the plot?
308
309    # FUNCTION:
310
311    # Settings:
312    # if (SPLUSLIKE) splusLikePlot(TRUE)
313
314    # Settings:
315    if (!is.timeSeries(x)) x = as.timeSeries(x)
316    Units = colnames(x)
317
318    # Labels:
319    if (labels) {
320        main = "Normal QQ Plot"
321        xlab = "Theoretical Quantiles"
322        ylab = "Sample Quantiles"
323    } else {
324        main = xlab = ylab = ""
325    }
326
327    X = x
328    DIM = dim(X)[2]
329
330    for (i in 1:DIM) {
331
332        # Get Data:
333        x = as.vector(as.matrix(X)[, i])
334        if (labels) main = Units[i]
335
336        # Standardized qqnorm():
337        y = (x-mean(x)) / sqrt(var(x))
338
339        # Further Settings:
340        y[abs(y) < span]
341        lim = c(-span, span)
342
343        # Plot qqnorm:
344        qqnorm(y, main = main, xlab = xlab, ylab = ylab,
345            xlim = lim, ylim = lim, col = col, ...)
346
347        # Add Line:
348        qqline(y, ...)
349
350        # Grid:
351        if (labels) grid()
352
353    }
354
355    # Return Value:
356    invisible(x)
357}
358
359
360# ------------------------------------------------------------------------------
361
362
363.ccfPlot <-
364function(x, y, lag.max = max(2, floor(10*log10(length(x)))),
365    type = c("correlation", "covariance", "partial"), labels = TRUE, ...)
366{
367    # A function implemented by Diethelm Wuertz
368
369    # Description:
370    #   Cross correlation function plot
371
372    # Arguments:
373    #   x - an univariate 'timeSeries' object
374    #   labels - a logical flag, by default true. Should a default
375    #       main title and labels addet to the plot?
376
377    # FUNCTION:
378
379    # Convert Type:
380    if (class(x) == "timeSeries") stopifnot(isUnivariate(x))
381    if (class(y) == "timeSeries") stopifnot(isUnivariate(y))
382    x = as.vector(x)
383    y = as.vector(y)
384
385    # Labels:
386    if (labels) {
387        main = "Crosscorrelation Function"
388        xlab = "lag"
389        ylab = "CCF"
390    } else {
391        main = xlab = ylab = ""
392    }
393
394    # Result:
395    # A copy from R's ccf - So you can use it also under SPlus:
396    X = cbind(x = x, y = y)
397    acf.out =  acf(X, lag.max = lag.max, plot = FALSE, type = type[1])
398    lag = c(rev(acf.out$lag[-1, 2, 1]), acf.out$lag[, 1, 2])
399    y = c(rev(acf.out$acf[-1, 2, 1]), acf.out$acf[, 1, 2])
400    acf.out$acf = array(y, dim = c(length(y), 1, 1))
401    acf.out$lag = array(lag, dim = c(length(y), 1, 1))
402    acf.out$snames = paste(acf.out$snames, collapse = " & ")
403    plot(acf.out, main = main, xlab = xlab, ylab = ylab, ...)
404
405    # Return Value:
406    invisible(acf.out)
407}
408
409
410################################################################################
411
412
413.stylizedFactsGUI <-
414function(x, mfrow = c(3, 3))
415{
416    # A function implemented by Diethelm Wuertz
417
418    # Description:
419    #   Opens a GUI for stylized facts
420
421    # FUNCTION:
422
423    # Refresh Code:
424    stylizedFactsRefreshCode <-
425function(...)
426    {
427        # Settings:
428        selectedAsset  = .tdSliderMenu(no = 1)
429        type = as.integer(.tdSliderMenu(obj.name = "stylizedFactsType"))
430        Unit = colnames(x)
431
432        # ACF Plot:
433        if (type == 1) {
434            if (selectedAsset == 0) {
435                par(mfrow = mfrow)
436                acfPlot(x)
437            } else {
438                par(mfrow = c(1, 1))
439                acfPlot(x[, selectedAsset])
440            }
441        }
442
443        # PACF Plot:
444        if (type == 2) {
445            if (selectedAsset == 0) {
446                par(mfrow = mfrow)
447                pacfPlot(x)
448            } else {
449                par(mfrow = c(1, 1))
450                pacfPlot(x[, selectedAsset])
451            }
452        }
453
454        # Volatility ACF Plot:
455        if (type == 3) {
456            if (selectedAsset == 0) {
457                par(mfrow = mfrow)
458                acfPlot(abs(x))
459            } else {
460                par(mfrow = c(1, 1))
461                acfPlot(abs(x[, selectedAsset]))
462            }
463        }
464
465        # Taylor Effect Plot:
466        if (type == 4) {
467            if (selectedAsset == 0) {
468                par(mfrow = mfrow)
469                teffectPlot(x)
470            } else {
471                par(mfrow = c(1, 1))
472                teffectPlot(x[, selectedAsset])
473            }
474        }
475
476        # Long Memory ACF:
477        if (type == 5) {
478            if (selectedAsset == 0) {
479                par(mfrow = mfrow)
480                .lmacfPlot(abs(x))
481            } else {
482                par(mfrow = c(1, 1))
483                .lmacfPlot(abs(x[, selectedAsset]))
484            }
485        }
486
487        # Lagged Autocorrelation Plot:
488        if (type == 6) {
489            if (selectedAsset == 0) {
490                par(mfrow = mfrow)
491                lacfPlot(x)
492            } else {
493                par(mfrow = c(1, 1))
494                lacfPlot(x[, selectedAsset])
495            }
496        }
497
498        # PDF plot on lin-log Scale:
499        if (type == 7) {
500            if (selectedAsset == 0) {
501                par(mfrow = mfrow)
502                .logpdfPlot(x)
503            } else {
504                par(mfrow = c(1, 1))
505                .logpdfPlot(x[, selectedAsset])
506            }
507        }
508
509        # PDF plot on log-log Scale:
510        if (type == 8) {
511            if (selectedAsset == 0) {
512                par(mfrow = mfrow)
513                .logpdfPlot(x, type = "log-log")
514            } else {
515                par(mfrow = c(1, 1))
516                .logpdfPlot(x[, selectedAsset], type = "log-log")
517            }
518        }
519
520        # Simple Normal Quantile-Quantile Plot
521        if (type == 9) {
522            if (selectedAsset == 0) {
523                par(mfrow = mfrow)
524                .qqgaussPlot(x, pch = 19)
525            } else {
526                par(mfrow = c(1, 1))
527                .qqgaussPlot(x[, selectedAsset], pch = 19)
528            }
529        }
530
531        # Scaling Law Plot:
532        if (type == 10) {
533            if (selectedAsset == 0) {
534                par(mfrow = mfrow)
535                scalinglawPlot(x, pch = 19)
536            } else {
537                par(mfrow = c(1, 1))
538                scalinglawPlot(x[, selectedAsset], pch = 19)
539            }
540        }
541
542    }
543
544    nAssets = dim(x)[2]
545
546    .tdSliderMenu(
547        stylizedFactsRefreshCode,
548
549        names       = c("Selected Asset"),
550        minima      = c(      0),
551        maxima      = c(      nAssets),
552        resolutions = c(      1),
553        starts      = c(      0),
554
555        but.functions = list(
556        function(...){
557                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "1")
558                stylizedFactsRefreshCode()},
559        function(...){
560                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "2")
561                stylizedFactsRefreshCode()},
562        function(...){
563                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "3")
564                stylizedFactsRefreshCode()},
565        function(...){
566                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "4")
567                stylizedFactsRefreshCode()},
568        function(...){
569                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "5")
570                stylizedFactsRefreshCode()},
571        function(...){
572                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "6")
573                stylizedFactsRefreshCode()},
574        function(...){
575                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "7")
576                stylizedFactsRefreshCode()},
577        function(...){
578                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "8")
579                stylizedFactsRefreshCode()},
580        function(...){
581                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "9")
582                stylizedFactsRefreshCode()},
583        function(...){
584                .tdSliderMenu(obj.name = "stylizedFactsType", obj.value = "10")
585                stylizedFactsRefreshCode()}
586        ),
587
588        but.names = c(
589            "1 ACF Function of Returns",
590            "2 Partial ACF of Returns",
591            "3 ACF of absolute Returns",
592            "4 Taylor Effect",
593            "5 Long Memory ACF of abs Returns",
594            "6 Lagged Autocorrelations",
595            "7 Lin-Log Tail Density Plot",
596            "8 Log-Log Lower Tail Density",
597            "9 Simple Normal QQ Plot",
598            "10 Scaling Law Plot"),
599
600        title = "Stylized Facts GUI"
601        )
602
603   .tdSliderMenu(obj.name = "type", obj.value = "1", no = 1)
604
605   # return Value:
606   invisible()
607}
608
609
610################################################################################
611
612