<!DOCTYPE html>
<!-- Generated by pkgdown: do not edit by hand --><html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
<meta charset="utf-8">
<meta http-equiv="X-UA-Compatible" content="IE=edge">
<meta name="viewport" content="width=device-width, initial-scale=1.0">
<title>Alternative climatologies and baselines • MHWdetection</title>
<!-- jquery --><script src="https://code.jquery.com/jquery-3.1.0.min.js" integrity="sha384-nrOSfDHtoPMzJHjVTdCopGqIqeYETSXhZDFyniQ8ZHcVy08QesyHcnOUpMpqnmWq" crossorigin="anonymous"></script><!-- Bootstrap --><link href="https://maxcdn.bootstrapcdn.com/bootstrap/3.3.7/css/bootstrap.min.css" rel="stylesheet" integrity="sha384-BVYiiSIFeK1dGmJRAkycuHAHRg32OmUcww7on3RYdg4Va+PmSTsz/K68vbdEjh4u" crossorigin="anonymous">
<script src="https://maxcdn.bootstrapcdn.com/bootstrap/3.3.7/js/bootstrap.min.js" integrity="sha384-Tc5IQib027qvyjSMfHjOMaLkfuWVxZxUPnCJA7l2mCWNIpG9mGCD8wGNIcPD7Txa" crossorigin="anonymous"></script><!-- Font Awesome icons --><link href="https://maxcdn.bootstrapcdn.com/font-awesome/4.6.3/css/font-awesome.min.css" rel="stylesheet" integrity="sha384-T8Gy5hrqNKT+hzMclPo118YTQO6cYprQmhrYwIiQ/3axmI1hQomh7Ud2hPOy8SP1" crossorigin="anonymous">
<!-- clipboard.js --><script src="https://cdnjs.cloudflare.com/ajax/libs/clipboard.js/1.7.1/clipboard.min.js" integrity="sha384-cV+rhyOuRHc9Ub/91rihWcGmMmCXDeksTtCihMupQHSsi8GIIRDG0ThDc3HGQFJ3" crossorigin="anonymous"></script><!-- pkgdown --><link href="../pkgdown.css" rel="stylesheet">
<script src="../jquery.sticky-kit.min.js"></script><script src="../pkgdown.js"></script><meta property="og:title" content="Alternative climatologies and baselines">
<meta property="og:description" content="This vignettes is about alternative climatologies and baselines.">
<meta name="twitter:card" content="summary">
<!-- mathjax --><script src="https://mathjax.rstudio.com/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script><!--[if lt IE 9]>
<script src="https://oss.maxcdn.com/html5shiv/3.7.3/html5shiv.min.js"></script>
<script src="https://oss.maxcdn.com/respond/1.4.2/respond.min.js"></script>
<![endif]-->
</head>
<body>
<div class="container template-article">
<header><div class="navbar navbar-default navbar-fixed-top" role="navigation">
<div class="container">
<div class="navbar-header">
<button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar">
<span class="icon-bar"></span>
<span class="icon-bar"></span>
<span class="icon-bar"></span>
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">MHWdetection</a>
<span class="label label-default" data-toggle="tooltip" data-placement="bottom" title="Released package">0.0.1</span>
</span>
</div>
<div id="navbar" class="navbar-collapse collapse">
<ul class="nav navbar-nav">
<li>
<a href="../index.html">
<span class="fa fa-home fa-lg"></span>
</a>
</li>
<li>
<a href="../reference/index.html">Functions</a>
</li>
<li class="dropdown">
<a href="#" class="dropdown-toggle" data-toggle="dropdown" role="button" aria-expanded="false">
Vignettes
<span class="caret"></span>
</a>
<ul class="dropdown-menu" role="menu">
<li>
<a href="../articles/Climatologies_and_baselines.html">Alternative climatologies and baselines</a>
</li>
<li>
<a href="../articles/Short_climatologies.html">How to create climatologies from short time series</a>
</li>
<li>
<a href="../articles/r_vs_python.html">Default MHW outputs in R and Python</a>
</li>
<li>
<a href="../articles/r_vs_python_additional.html">Additional functionality for MHW outputs in R and Python</a>
</li>
<li>
<a href="../articles/r_vs_python_arguments.html">Arguments for MHW outputs in R and Python</a>
</li>
</ul>
</li>
<li>
<a href="../news/index.html">Changelog</a>
</li>
</ul>
<ul class="nav navbar-nav navbar-right"></ul>
</div>
<!--/.nav-collapse -->
</div>
<!--/.container -->
</div>
<!--/.navbar -->
</header><div class="row">
<div class="col-md-9 contents">
<div class="page-header toc-ignore">
<h1>Alternative climatologies and baselines</h1>
<h4 class="author">AJ Smit</h4>
<h4 class="date">2018-05-24</h4>
<div class="hidden name"><code>Climatologies_and_baselines.Rmd</code></div>
</div>
<div id="background" class="section level1">
<h1 class="hasAnchor">
<a href="#background" class="anchor"></a>Background</h1>
<p>A time series of temperatures taken at a specific locality is a sample of the thermal environment. Because it is a sample, which is used to approximate the true state of the environment, the representivity of the data is subject to certain limitations that have an influence on their accuracy and precision—in other words, on how well the data represent reality. Here we will not concern ourselves with matters of precision and accuracy; rather, we shall deal with ways to produce climatologies from a time series of daily thermal measurements.</p>
<p>A climatological mean is a multiyear average over a base period of typically 30 years duration. Sometimes the climatology is also called a ‘climate normal.’ It predicts the temperatures over an annual cycle that will most likely be experienced at a particlar place and time.</p>
<p>Another concept that we use here instead of a climatology is a baseline… <definition to follow>.</definition></p>
<div class="figure">
<img src="fig/heatwaveR_v3.svg" alt="Figure 1. Processing pipeline in the detect() algorithm with emphasis on the climatology and baseline procedures. The dark black arrows indicate the options that are currently available." style="width:100.0%"><p class="caption"><strong>Figure 1.</strong> Processing pipeline in the <code>detect()</code> algorithm with emphasis on the climatology and baseline procedures. The dark black arrows indicate the options that are currently available.</p>
</div>
<div id="the-default-climatology-and-threshold" class="section level2">
<h2 class="hasAnchor">
<a href="#the-default-climatology-and-threshold" class="anchor"></a>The default climatology and threshold</h2>
<p>Central to the functioning of the <strong>heatwaveR</strong> detect algogorithm is the climatology, which is the reference against which the extremes are calculated. For historical reasons—which also turn out to make good practical sense—it has become standard practice to base a climatology on 30 years’ worth of data, updated at the start of each new decade (the most recent being 1981-2010). For the <strong>heatwaveR</strong> package we recommend the same. Because heat wave detection is based on daily data, a daily climatology is required. That is, the mean value for each of the days in a 365 (non-leap) or 366 (leap) day long year. In the case of gridded data, we also require that the climatology is at the same spatial resolution as the SST time series within which the extreme events are sought. In the case of the dOISST v2. data, it is independently calculated at each of the 1/4° grid cells, and for station data it must be representative of that station and preferably also obtained with the same instrument that yielded the seawater temperature time series entering the event detection algorithm (in fact, the algorithm will by default use the daily time series to derive the climatology).</p>
<p>In a daily climatology, there will be one mean temperature for each day of the year (365 or 366 days, depending on how leap years are dealt with). How the daily mean is calculated can vary somewhat from application to application. The <code>detect()</code> function applies a window centered on a sliding Julian day for the pooling of values that will form the basis for calculation of the climatology and threshold percentile. We use a default of 5 days, which gives a window width of 11 days centered on the 6th day of the series of 11 days. For example, the mean temperature for Julian day 6 will be comprised of all temperatures measured from Julian days 1 to 11, across 30 years, i.e. 330 temperature measurements are summarised in this climatological daily mean. A sliding window centered on Julian day 7 will comprise all temperatures measured on Julian days 2 to 12 over the 30 years; etc. Currently we do not provide the option to calculate a weighted mean; rather, each day receives the same weighting as any other day, irrespective of how far the day is from the one on which the sliding window is centered. Climatologies produced with sliding windows of width 1, 11, 31, and 61 days are presented in Figure 1. In the figure we see that the wider the sliding window becomes the less day-to-day noise is retained in the climatology.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">library</span>(tidyverse) <span class="co"># misc. data processing conveniences</span>
<span class="kw">library</span>(lubridate) <span class="co"># for working with dates</span>
<span class="kw">library</span>(heatwaveR)</code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># Using a built-in data set</span>
<span class="co"># smoothPercentile = FALSE</span>
res.<span class="fl">1.</span>F <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>),
<span class="dt">windowHalfWidth =</span> <span class="dv">0</span>, <span class="dt">smoothPercentile =</span> F, <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>)
res.<span class="fl">11.</span>F <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>),
<span class="dt">smoothPercentile =</span> F, <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>)
res.<span class="fl">31.</span>F <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>),
<span class="dt">windowHalfWidth =</span> <span class="dv">15</span>, <span class="dt">smoothPercentile =</span> F, <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>)
res.<span class="fl">61.</span>F <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>),
<span class="dt">windowHalfWidth =</span> <span class="dv">30</span>, <span class="dt">smoothPercentile =</span> F, <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>)</code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">res.<span class="fl">1.</span>F <-<span class="st"> </span>res.<span class="fl">1.</span>F <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"1"</span>)
res.<span class="fl">11.</span>F <-<span class="st"> </span>res.<span class="fl">11.</span>F <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"11"</span>)
res.<span class="fl">31.</span>F <-<span class="st"> </span>res.<span class="fl">31.</span>F <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"31"</span>)
res.<span class="fl">61.</span>F <-<span class="st"> </span>res.<span class="fl">61.</span>F <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"61"</span>)
no_smooth <-<span class="st"> </span><span class="kw">rbind</span>(res.<span class="fl">1.</span>F, res.<span class="fl">11.</span>F, res.<span class="fl">31.</span>F, res.<span class="fl">61.</span>F)
plt1 <-<span class="st"> </span><span class="kw">ggplot</span>(<span class="dt">data =</span> no_smooth,
<span class="kw">aes</span>(<span class="dt">x =</span> doy, <span class="dt">y =</span> seas)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="kw">aes</span>(<span class="dt">colour =</span> width, <span class="dt">size =</span> width)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="kw">aes</span>(<span class="dt">y =</span> thresh, <span class="dt">colour =</span> width, <span class="dt">size =</span> width)) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_colour_manual</span>(<span class="dt">name =</span> <span class="st">"Width"</span>,
<span class="dt">values =</span> <span class="kw">c</span>(<span class="st">"grey50"</span>, <span class="st">"orange"</span>, <span class="st">"purple"</span>, <span class="st">"black"</span>)) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_size_manual</span>(<span class="dt">values =</span> <span class="kw">c</span>(<span class="fl">0.4</span>, <span class="fl">0.4</span>, <span class="fl">0.4</span>, <span class="fl">0.8</span>), <span class="dt">guide =</span> <span class="ot">FALSE</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">labs</span>(<span class="dt">x =</span> <span class="st">"Julian day"</span>, <span class="dt">y =</span> <span class="st">"SST (°C)"</span>, <span class="dt">title =</span> <span class="st">"Annual cycle climatology and threshold"</span>,
<span class="dt">subtitle =</span> <span class="st">"1983-01-01 to 2012-12-31 reference period; smooth off"</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">theme_publ</span>()
<span class="kw">ggsave</span>(plt1, <span class="dt">file =</span> <span class="st">"climat_1.png"</span>, <span class="dt">width =</span> <span class="dv">6</span>, <span class="dt">height =</span> <span class="fl">2.8</span>, <span class="dt">dpi =</span> <span class="dv">120</span>)</code></pre></div>
<div class="figure">
<img src="fig/climat_1.png" alt="Figure 2. Daily climatologies of the mean (lower lines) and 90th percentile (top lines) produced by the detect() function, with various windowHalfWidth values and with smoothPercentile switched off." style="width:80.0%"><p class="caption"><strong>Figure 2.</strong> Daily climatologies of the mean (lower lines) and 90th percentile (top lines) produced by the <code>detect()</code> function, with various <code>windowHalfWidth</code> values and with <code>smoothPercentile</code> switched off.</p>
</div>
<p>To produce smoother climatologies, we can additionally enable a moving average smoother. The result is shown in Figure 3.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co">#smoothPercentile = TRUE</span>
res.<span class="fl">11.</span>T <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>)
res.<span class="fl">31.</span>T <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>, <span class="dt">windowHalfWidth =</span> <span class="dv">15</span>)
res.<span class="fl">61.</span>T <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">TRUE</span>, <span class="dt">windowHalfWidth =</span> <span class="dv">30</span>)</code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">res.<span class="fl">11.</span>T <-<span class="st"> </span>res.<span class="fl">11.</span>T <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"11"</span>)
res.<span class="fl">31.</span>T <-<span class="st"> </span>res.<span class="fl">31.</span>T <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"31"</span>)
res.<span class="fl">61.</span>T <-<span class="st"> </span>res.<span class="fl">61.</span>T <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">width =</span> <span class="st">"61"</span>)
smooth <-<span class="st"> </span><span class="kw">rbind</span>(res.<span class="fl">11.</span>T, res.<span class="fl">31.</span>T, res.<span class="fl">61.</span>T)
plt2 <-<span class="st"> </span><span class="kw">ggplot</span>(<span class="dt">data =</span> smooth,
<span class="kw">aes</span>(<span class="dt">x =</span> doy, <span class="dt">y =</span> seas)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="kw">aes</span>(<span class="dt">colour =</span> width, <span class="dt">size =</span> width)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="kw">aes</span>(<span class="dt">y =</span> thresh, <span class="dt">colour =</span> width, <span class="dt">size =</span> width)) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_colour_manual</span>(<span class="dt">name =</span> <span class="st">"Width"</span>,
<span class="dt">values =</span> <span class="kw">c</span>(<span class="st">"orange"</span>, <span class="st">"purple"</span>, <span class="st">"black"</span>)) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_size_manual</span>(<span class="dt">values =</span> <span class="kw">c</span>(<span class="fl">0.4</span>, <span class="fl">0.4</span>, <span class="fl">0.8</span>), <span class="dt">guide =</span> <span class="ot">FALSE</span>) <span class="op">+</span><span class="st"> </span>
<span class="st"> </span><span class="kw">labs</span>(<span class="dt">x =</span> <span class="st">"Julian day"</span>, <span class="dt">y =</span> <span class="st">"SST (°C)"</span>, <span class="dt">title =</span> <span class="st">"Annual cycle climatology and threshold"</span>,
<span class="dt">subtitle =</span> <span class="st">"1983-01-01 to 2012-12-31 reference period; smooth on"</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">theme_publ</span>()
<span class="kw">ggsave</span>(plt2, <span class="dt">file =</span> <span class="st">"climat_2.png"</span>, <span class="dt">width =</span> <span class="dv">6</span>, <span class="dt">height =</span> <span class="fl">2.8</span>, <span class="dt">dpi =</span> <span class="dv">120</span>)</code></pre></div>
<div class="figure">
<img src="fig/climat_2.png" alt="Figure 3. Daily climatologies of the mean (lower lines) and 90th percentile (upper lines) produced by the detect() function, with various windowHalfWidth values and with smoothPercentile at the default value of 31 days." style="width:80.0%"><p class="caption"><strong>Figure 3.</strong> Daily climatologies of the mean (lower lines) and 90th percentile (upper lines) produced by the <code>detect()</code> function, with various <code>windowHalfWidth</code> values and with <code>smoothPercentile</code> at the default value of 31 days.</p>
</div>
<div id="fourier-analysis-of-time-series-to-construct-a-smooth-climatology" class="section level3">
<h3 class="hasAnchor">
<a href="#fourier-analysis-of-time-series-to-construct-a-smooth-climatology" class="anchor"></a>Fourier analysis of time series to construct a smooth climatology</h3>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">library</span>(fda)
clm <-<span class="st"> </span>res.<span class="fl">1.</span>F <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span>as_tibble
b7 <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/fda/topics/create.fourier.basis">create.fourier.basis</a></span>(<span class="dt">rangeval =</span> <span class="kw">range</span>(clm<span class="op">$</span>doy), <span class="dt">nbasis =</span> <span class="dv">7</span>)
<span class="co"># plot(b7)</span>
<span class="co"># b7.val <- as_tibble(eval.basis(clm$doy, basisobj = b7))</span>
<span class="co"># b7.val <- cbind(doy = clm$doy, b7.val) %>% </span>
<span class="co"># gather(key = "basis", value = "value", -doy) %>% </span>
<span class="co"># as_tibble()</span>
<span class="co"># create smooth seasonal climatology</span>
b7.seas.smth <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/fda/topics/smooth.basis">smooth.basis</a></span>(<span class="dt">argvals =</span> clm<span class="op">$</span>doy, <span class="dt">y =</span> clm<span class="op">$</span>seas, <span class="dt">fdParobj =</span> b7)<span class="op">$</span>fd
clm.seas.smth <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/fda/topics/eval.fd">eval.fd</a></span>(clm<span class="op">$</span>doy, b7.seas.smth)
<span class="co"># create smooth threshold</span>
b7.thresh.smth <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/fda/topics/smooth.basis">smooth.basis</a></span>(<span class="dt">argvals =</span> clm<span class="op">$</span>doy, <span class="dt">y =</span> clm<span class="op">$</span>thresh, <span class="dt">fdParobj =</span> b7)<span class="op">$</span>fd
clm.thresh.smth <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/fda/topics/eval.fd">eval.fd</a></span>(clm<span class="op">$</span>doy, b7.thresh.smth)
temps <-<span class="st"> </span><span class="kw">as_tibble</span>(<span class="kw">data.frame</span>(<span class="dt">doy =</span> clm<span class="op">$</span>doy,
<span class="dt">seas.unsmth_1 =</span> clm<span class="op">$</span>seas,
<span class="dt">seas.smth_11 =</span> res.<span class="fl">11.</span>T<span class="op">$</span>seas,
<span class="dt">clim.fourier_7 =</span> clm.seas.smth,
<span class="dt">thresh.unsmth_1 =</span> clm<span class="op">$</span>thresh,
<span class="dt">thresh.smth_11 =</span> res.<span class="fl">11.</span>T<span class="op">$</span>thresh,
<span class="dt">thresh.fourier_7 =</span> clm.thresh.smth))</code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">library</span>(forcats)
plt3 <-<span class="st"> </span>temps <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">gather</span>(<span class="dt">key =</span> <span class="st">"type"</span>, <span class="dt">value =</span> <span class="st">"temperature"</span>, <span class="op">-</span>doy) <span class="op">%>%</span><span class="st"> </span>
<span class="st"> </span><span class="kw">ggplot</span>(<span class="kw">aes</span>(<span class="dt">x =</span> doy, <span class="dt">y =</span> temperature)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="kw">aes</span>(<span class="dt">colour =</span> <span class="kw"><a href="http://www.rdocumentation.org/packages/forcats/topics/fct_inorder">fct_inorder</a></span>(type), <span class="dt">size =</span> <span class="kw"><a href="http://www.rdocumentation.org/packages/forcats/topics/fct_inorder">fct_inorder</a></span>(type))) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_colour_manual</span>(<span class="dt">name =</span> <span class="st">"Smooth"</span>,
<span class="dt">values =</span> <span class="kw">c</span>(<span class="st">"grey50"</span>, <span class="st">"purple"</span>, <span class="st">"black"</span>, <span class="st">"grey50"</span>, <span class="st">"purple"</span>, <span class="st">"black"</span>),
<span class="dt">labels =</span> <span class="kw">c</span>(<span class="st">"1, unsmooth"</span>, <span class="st">"11, smooth"</span>, <span class="st">"Fourier"</span>),
<span class="dt">breaks =</span> <span class="kw">c</span>(<span class="st">"seas.unsmth_1"</span>, <span class="st">"seas.smth_11"</span>, <span class="st">"clim.fourier_7"</span>)) <span class="op">+</span>
<span class="st"> </span><span class="kw">scale_size_manual</span>(<span class="dt">values =</span> <span class="kw">c</span>(<span class="fl">0.4</span>, <span class="fl">0.4</span>, <span class="fl">0.8</span>, <span class="fl">0.4</span>, <span class="fl">0.4</span>, <span class="fl">0.8</span>), <span class="dt">guide =</span> <span class="ot">FALSE</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">labs</span>(<span class="dt">x =</span> <span class="st">"Julian day"</span>, <span class="dt">y =</span> <span class="st">"SST (°C)"</span>, <span class="dt">title =</span> <span class="st">"Annual cycle climatology and threshold"</span>,
<span class="dt">subtitle =</span> <span class="st">"1983-01-01 to 2012-12-31 reference period"</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">theme_publ</span>()
<span class="kw">ggsave</span>(plt3, <span class="dt">file =</span> <span class="st">"climat_3.png"</span>, <span class="dt">width =</span> <span class="dv">6</span>, <span class="dt">height =</span> <span class="fl">2.8</span>, <span class="dt">dpi =</span> <span class="dv">120</span>)
<span class="co"># plot(b7.seas.smth[1:7])</span>
<span class="co"># b7.PCA = pca.fd(b7.seas.smth, nharm=6)</span></code></pre></div>
<div class="figure">
<img src="fig/climat_3.png" alt="Figure 4. Daily climatologies of the mean (lower lines) and 90th percentile (upper lines) produced by the detect() function; 1, unsmooth indicates that a windowHalfWidth of 0 days had been set, and that the data had not been subjected to smoothPercentile; 11, smooth is the default settings of the detect() function; and Fourier indicates that an alternative climatology had been produced through the application of a Fourier analysis with 7 basis points." style="width:80.0%"><p class="caption"><strong>Figure 4.</strong> Daily climatologies of the mean (lower lines) and 90th percentile (upper lines) produced by the <code>detect()</code> function; “1, unsmooth” indicates that a <code>windowHalfWidth</code> of 0 days had been set, and that the data had not been subjected to <code>smoothPercentile</code>; “11, smooth” is the default settings of the <code>detect()</code> function; and “Fourier” indicates that an alternative climatology had been produced through the application of a Fourier analysis with 7 basis points.</p>
</div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># library(TSA)</span>
<span class="co"># # compute the Fourier Transform</span>
<span class="co"># p <- periodogram(y)</span>
<span class="co"># dd <- data.frame(freq = p$freq, spec = p$spec)</span>
<span class="co"># order <- dd[order(-dd$spec), ]</span>
<span class="co"># top6 <- head(order, 6)</span>
<span class="co"># # display the 6 highest "power" frequencies</span>
<span class="co"># top6</span>
<span class="co"># # convert frequency to time periods</span>
<span class="co"># time <- 1/top6$f</span>
<span class="co"># time</span></code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># fft(y)</span>
<span class="co"># spectrum(y)</span></code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># require(graphics)</span>
<span class="co"># spec.y <- spec.pgram(y)</span>
<span class="co"># plot(spec.y, plot.type = "coherency")</span>
<span class="co"># plot(spec.y, plot.type = "phase")</span>
<span class="co"># library(GeneCycle)</span>
<span class="co"># GeneCycle::periodogram(y)</span></code></pre></div>
</div>
</div>
<div id="using-a-custom-baseline" class="section level2">
<h2 class="hasAnchor">
<a href="#using-a-custom-baseline" class="anchor"></a>Using a custom baseline</h2>
<p>Should one want to use a different baseline to calculate the climatology against which the events should be detected, this is now possible.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># Again, using the built-in time series</span>
<span class="co"># The Hobday et al. (2016) default definition applied above is available in</span>
<span class="co"># `res.11.T`; we repeat it her and return the original time series with the</span>
<span class="co"># climatology so as to have data in the specified format as required by the</span>
<span class="co"># detect function, `detect_event()`:</span>
res.<span class="fl">11.</span>T.full <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(sst_WA, <span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">FALSE</span>)
<span class="co"># Calculate an anomaly:</span>
daily.dat <-<span class="st"> </span>res.<span class="fl">11.</span>T.full <span class="op">%>%</span>
<span class="st"> </span><span class="kw">mutate</span>(<span class="dt">anom =</span> temp <span class="op">-</span><span class="st"> </span><span class="kw">mean</span>(temp, <span class="dt">na.rm =</span> <span class="ot">TRUE</span>),
<span class="dt">num =</span> <span class="kw">seq</span>(<span class="dv">1</span><span class="op">:</span><span class="kw">length</span>(temp)))
<span class="co"># Let us remove the long-term non-linear trend from the WA anomaly data using a </span>
<span class="co"># generalised additive model; first, define the GAM:</span>
<span class="kw">library</span>(broom)
mod.gam <-<span class="st"> </span><span class="kw">gam</span>(temp <span class="op">~</span><span class="st"> </span><span class="kw">s</span>(num, <span class="dt">k =</span> <span class="dv">20</span>), <span class="dt">data =</span> daily.dat)
daily.dat2 <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/broom/topics/augment">augment</a></span>(mod.gam)
daily.dat2 <-<span class="st"> </span><span class="kw">cbind</span>(daily.dat2, daily.dat[, <span class="kw">c</span>(<span class="st">"doy"</span>, <span class="st">"t"</span>, <span class="st">"anom"</span>)])
raw.clim <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(daily.dat2, <span class="dt">x =</span> t, <span class="dt">y =</span> anom,
<span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">FALSE</span>)
detrended.clim <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/ts2clm">ts2clm</a></span>(daily.dat2, <span class="dt">x =</span> t, <span class="dt">y =</span> .resid,
<span class="dt">climatologyPeriod =</span> <span class="kw">c</span>(<span class="st">"1983-01-01"</span>, <span class="st">"2012-12-31"</span>), <span class="dt">clmOnly =</span> <span class="ot">FALSE</span>)
<span class="co"># Plot the two climatologies (raw vs. detrended):</span>
plt4 <-<span class="st"> </span><span class="kw">ggplot</span>(<span class="dt">data =</span> <span class="kw">filter</span>(raw.clim, t <span class="op"><</span><span class="st"> "1983-01-01"</span>), <span class="kw">aes</span>(<span class="dt">x =</span> doy, <span class="dt">y =</span> seas)) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="dt">colour =</span> <span class="st">"black"</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">geom_line</span>(<span class="dt">data =</span> <span class="kw">filter</span>(detrended.clim, t <span class="op"><</span><span class="st"> "1983-01-01"</span>),
<span class="dt">colour =</span> <span class="st">"red"</span>, <span class="dt">alpha =</span> <span class="fl">0.6</span>) <span class="op">+</span>
<span class="st"> </span><span class="kw">theme_publ</span>()
<span class="kw">ggsave</span>(plt4, <span class="dt">file =</span> <span class="st">"climat_4.png"</span>, <span class="dt">width =</span> <span class="dv">12</span>, <span class="dt">height =</span> <span class="dv">6</span>, <span class="dt">dpi =</span> <span class="dv">120</span>)
<span class="co"># ... there is a tiny, almost imperceptable difference between the two...</span>
<span class="co"># Calculate events in the original WA SST (anomaly) relative to a climatology </span>
<span class="co"># created from the detrended WA SST data:</span>
raw.res <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/detect_event">detect_event</a></span>(raw.clim, <span class="dt">x =</span> t, <span class="dt">y =</span> anom,
<span class="dt">seasClim =</span> seas,
<span class="dt">threshClim =</span> thresh)
detrended.res <-<span class="st"> </span><span class="kw"><a href="http://www.rdocumentation.org/packages/heatwaveR/topics/detect_event">detect_event</a></span>(raw.clim, <span class="dt">x =</span> t, <span class="dt">y =</span> anom,
<span class="dt">seasClim =</span> detrended.clim<span class="op">$</span>seas,
<span class="dt">threshClim =</span> detrended.clim<span class="op">$</span>thresh)</code></pre></div>
<div class="figure">
<img src="fig/climat_4.png" alt="Figure 5. Daily climatologies of the WA SST (anomalies). The black line is the climatology prepared from the raw data, and the red line represents the daily SST anomalies that had been detrended using a generalised additive model. The lines plot practically on top of each other, so differences are barely perceptable; yet, despite these apparently small differences, there are large consequences for the events that are detected, as can be seen in the printout of the top 10 events, as seen below." style="width:80.0%"><p class="caption"><strong>Figure 5.</strong> Daily climatologies of the WA SST (anomalies). The black line is the climatology prepared from the raw data, and the red line represents the daily SST anomalies that had been detrended using a generalised additive model. The lines plot practically on top of each other, so differences are barely perceptable; yet, despite these apparently small differences, there are large consequences for the events that are detected, as can be seen in the printout of the top 10 events, as seen below.</p>
</div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">raw.res</code></pre></div>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">detrended.res</code></pre></div>
<p>The events that get detected are somewhat different now. Notice that in <code>detrended.clim</code>, <code>seas</code> and <code>thresh</code> are unfortunately the ones that were calculted together with raw (anomaly) time series, not the <code>seas</code> and <code>thresh</code> that were actually used and included against climSeas and climThresh. We will fix that.</p>
</div>
<div id="using-a-custom-climatology" class="section level2">
<h2 class="hasAnchor">
<a href="#using-a-custom-climatology" class="anchor"></a>Using a custom climatology</h2>
<p>If one would rather supply a pre-calculated climatology, and not calculate one from a baseline, this is also now possible.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># Pull out the climatologies from the built in time series</span>
<span class="co"># ts1_clim <- detect_evets2clm(sst_WA, climatologyPeriod = c("1983-01-01", "2012-12-31"))$clim %>%</span>
<span class="co"># dplyr::select(doy, seas, thresh, var_clim_year) %>%</span>
<span class="co"># base::unique() %>%</span>
<span class="co"># dplyr::arrange(doy)</span>
<span class="co"># ts2_clim <- detect_evets2clm(sst_WA, climatologyPeriod = c("1983-01-01", "2012-12-31"))$clim %>%</span>
<span class="co"># dplyr::select(doy, seas, thresh, var_clim_year) %>%</span>
<span class="co"># base::unique() %>%</span>
<span class="co"># dplyr::arrange(doy)</span>
<span class="co"># ts3_clim <- detect_evets2clm(sst_WA, climatologyPeriod = c("1983-01-01", "2012-12-31"))$clim %>%</span>
<span class="co"># dplyr::select(doy, seas, thresh, var_clim_year) %>%</span>
<span class="co"># base::unique() %>%</span>
<span class="co"># dplyr::arrange(doy)</span>
<span class="co"># </span>
<span class="co"># # Calculate events in one time series with a climatology from another</span>
<span class="co"># res_clim_diff <- detect_event(ts1, alt_clim = TRUE, alt_clim_data = ts2_clim)</span></code></pre></div>
<p>Note how the climatologies are the same below:</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># head(ts2_clim)</span>
<span class="co"># head(res_clim_diff$clim)</span></code></pre></div>
<p>Should the provided climatology be to low, the function will return a warning letting the user know this.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="co"># # Calculate events in one time series with a climatology from another</span>
<span class="co"># res_clim_diff <- detect_event(ts1, alt_clim = TRUE, alt_clim_data = ts3_clim)</span></code></pre></div>
</div>
</div>
</div>
<div class="col-md-3 hidden-xs hidden-sm" id="sidebar">
<div id="tocnav">
<h2 class="hasAnchor">
<a href="#tocnav" class="anchor"></a>Contents</h2>
<ul class="nav nav-pills nav-stacked">
<li>
<a href="#background">Background</a><ul class="nav nav-pills nav-stacked">
<li><a href="#the-default-climatology-and-threshold">The default climatology and threshold</a></li>
<li><a href="#using-a-custom-baseline">Using a custom baseline</a></li>
<li><a href="#using-a-custom-climatology">Using a custom climatology</a></li>
</ul>
</li>
</ul>
</div>
</div>
</div>
<footer><div class="copyright">
<p>Developed by Robert W. Schlegel, Eric C. J. Oliver, Alistair J. Hobday, Albertus J. Smit.</p>
</div>
<div class="pkgdown">
<p>Site built with <a href="http://pkgdown.r-lib.org/">pkgdown</a>.</p>
</div>
</footer>
</div>
</body>
</html>