TheilSen(mydata, pollutant = "nox", deseason = FALSE, type = "default",
avg.time = "month", statistic = "mean", percentile = NA, data.thresh = 0,
alpha = 0.05, dec.place = 2, xlab = "year", lab.frac = 0.99, lab.cex = 0.8,
x.relation = "same", y.relation = "same", data.col = "cornflowerblue",
line.col = "red", text.col = "darkgreen", cols = NULL, auto.text = TRUE,
autocor = FALSE, slope.percent = FALSE, date.breaks = 7,...) MannKendall(mydata,...)
date and at least one other parameter for which a
trend test is required; typically (but not necessarily) a
pollutant.TRUE the function stl is used (seasonal
trend decomposition using loess). Note that if
TRUE missing data are first linearly interpolated
because stl ctype determines how the data are split
i.e. conditioned, and then plotted. The default is will
produce a single plot using the entire data. Type can be
one of the built-in types as detailed in cutData
e.g. timeAverage for more
details.statistic = "percentile" is chosen.avg.time. A value
of zero means that all available data will be used in a
particular period regardless if of the number of values
available. Conversely, a value of 100"year"."greyscale".TRUE (default) or
FALSE. If TRUE titles and axis labels will
automatically try and format pollutant names and units
properly e.g. by subscripting the FALSE. Generally, accounting for autocorrelation
increases the uncertainty of the trend estimate ---
sometimes by a large amount.FALSE and the slope is
expressed as an average units/year change e.g. ppb.
Percentage changes can often be confusing and should cutData and lattice:xyplot. For example,
TheilSen passes the option hemisphere =
"southern" on to cutData to provide southern
(rather than defaTheilSen also
returns an object of class ``openair''. The object includes
three main components: call, the command used to
generate the plot; data, the data frame of
summarised information used to make the plot; and
plot, the plot itself. If retained, e.g. using
output <- TheilSen(mydata, "nox"), this output can
be used to recover the data, reproduce or rework the
original plot or undertake further analysis.An openair output can be manipulated using a number of
generic operations, including print, plot and
summary. See openair.generics for
further details.
The data component of the TheilSen output
includes two subsets: main.data, the monthly data
res2 the trend statistics. For output <-
TheilSen(mydata, "nox"), these can be extracted as
object$data$main.data and object$data$res2,
respectively.
Note: In the case of the intercept, it is assumed the y-axis crosses the x-axis on 1/1/1970.
TheilSen function provides a collection of
functions to analyse trends in air pollution data. The
Mann-Kendall test is a commonly used test in environmental
sciences to detect the presence of a trend. It is often
used with the Theil-Sen (or just Sen) estimate of slope.
See references. The TheilSen function is flexible
in the sense that it can be applied to data in many ways
e.g. by day of the week, hour of day and wind direction.
This flexibility makes it much easier to draw inferences
from data e.g. why is there a strong downward trend in
concentration from one wind sector and not another, or why
trends on one day of the week or a certain time of day are
unexpected.For data that are strongly seasonal, perhaps from a
background site, or a pollutant such as ozone, it will be
important to deseasonalise the data (using the option
deseason = TRUE.Similarly, for data that increase,
then decrease, or show sharp changes it may be better to
use smoothTrend.
Note! that since version 0.5-11 openair uses Theil-Sen to
derive the p values also. This is to ensure there is
consistency between the calculated p value and other trend
parameters i.e. slope estimates and uncertainties. This
change may slightly affect some of the p-estimates
previously given by openair because the p estimates are now
calculated using bootstrap resampling by default and
previously they were not. However, users can still for the
moment call the TheilSen function using
MannKendall.
Note that the symbols shown next to each trend estimate relate to how statistically significant the trend estimate is: p $<$ 0.001="***," p="" $<$="" 0.01="**," 0.05="*" and="" 0.1="$+$.
Some of the code used in TheilSen is based on that
from Rand Wilcox autocor = TRUE (Kunsch, 1989). We follow the
suggestion of Kunsch (1989) of setting the block length to
n(1/3) where n is the length of the time series.
The slope estimate and confidence intervals in the slope are plotted and numerical information presented.
Hirsch, R. M., Slack, J. R., Smith, R. A., 1982. Techniques of trend analysis for monthly water-quality data. Water Resources Research 18 (1), 107-121.
Kunsch, H. R., 1989. The jackknife and the bootstrap for general stationary observations. Annals of Statistics 17 (3), 1217-1241.
Sen, P. K., 1968. Estimates of regression coefficient based on Kendall's tau. Journal of the American Statistical Association 63(324).
Theil, H., 1950. A rank invariant method of linear and polynomial regression analysis, i, ii, iii. Proceedings of the Koninklijke Nederlandse Akademie Wetenschappen, Series A - Mathematical Sciences 53, 386-392, 521-525, 1397-1412.
...
smoothTrend for a flexible approach to
estimating trends using nonparametric regression. The
smoothTrend function is suitable for cases where
trends are not monotonic and is probably better for
exploring the shape of trends.# load example data from package
data(mydata)
# trend plot for nox
TheilSen(mydata, pollutant = "nox")
# trend plot for ozone with p=0.01 i.e. uncertainty in slope shown at
# 99 \% confidence interval
TheilSen(mydata, pollutant = "o3", ylab = "o3 (ppb)", alpha = 0.01)
# trend plot by each of 8 wind sectors
TheilSen(mydata, pollutant = "o3", type = "wd", ylab = "o3 (ppb)")
# and for a subset of data (from year 2000 onwards)
TheilSen(select.by.date(mydata, year = 2000:2005), pollutant = "o3", ylab = "o3 (ppb)")Run the code above in your browser using DataLab