---
output: 
  html_document:
    css: help.css

params:
  compartment: sediment
  detGroup: Organo-metals
---

```{r setup, include = FALSE, echo = FALSE}

library(knitr)
knit_engines$set(asis = function(options) {
  if (options$echo && options$eval) knit_child(text = options$code)
})

options(width = 110)

periodStart <- 2010
periodEnd <- 2015
assessmentYear <- 2017

source("report functions.R")

load(file.path("meta analysis objects", "regional data.RData"))
assessment <- get(params$compartment)
assessment <- assessment[[params$detGroup]]
purpose <- assessment$purpose

isAC <- !is.null(assessment$AC)

detAC <- switch(
  params$compartment,
  sediment = switch(
    params$detGroup,
    Metals = "Effects Range Low",
    Chlorobiphenyls = "Environmental Assessment Criterion",
    "PAH (parent)" = "Effects Range Low",
    NULL
  ),
  biota = switch(
    params$detGroup,
    Metals = "European Commission food standard",
    Chlorobiphenyls = "Environmental Assessment Criterion",
    "PAH (parent)" = "Environmental Assessment Criterion",
    NULL
  )
)

if (isAC) {
  detAClabel <- switch(
    detAC, 
    "Effects Range Low" = "ERL",
    "European Commission food standard" = "EC",
    "Environmental Assessment Criterion" = "EAC"
  )
}

detLabel <- switch(
  params$detGroup, 
  Metals = "metal", 
  Chlorobiphenyls = "CB",
  "PAH (parent)" = "PAH",
  "Organo-bromines" = "BDE",
  "Organo-metals" = "organo-metal"
)


detFixedEffect <- params$detGroup %in% "Metals"

if (params$compartment == "biota" && detFixedEffect) {
  greenID <- "orange" 
  isOrange <- FALSE
} else {
  greenID <- "green"
  isOrange <- any(assessment$status$colour %in% "orange")
}
  
isBlack <- any(assessment$status$colour %in% "black")

detUnit <- switch(
  params$compartment, 
  sediment = switch(
    params$detGroup, 
    Metals = "mg/kg dw", 
    Chlorobiphenyls = "ug/kg dw",
    "PAH (parent)" = "ug/kg dw",
    "Organo-bromines" = "ug/kg dw",
    "Organo-metals" = "ug/kg dw"
  ),
  biota = switch(
    params$detGroup,
    Metals = "ug/kg ww", 
    Chlorobiphenyls = "ug/kg ww",
    "PAH (parent)" = "ug/kg ww",
    "Organo-bromines" = "ug/kg ww",
    "Organo-metals" = "ug/kg ww"
  )
)


```

## MIME regional assessment of status and trends in `r detLabel` concentrations in `r params$compartment`

### Overview

<p class = "startbig">
The `r assessmentYear` MIME assessment describes the trends and status of contaminant concentrations in biota and sediment at monitoring stations in the OSPAR area. Assessments are made for a large number of time series, each of a single contaminant in a single species (for biota) at a single monitoring station.  This document is one of a series that synthesises the results of the individual time series to assess status and trends at the MIME regional level.  In particular, it considers `r detLabel` concentrations in `r params$compartment`, where the `r detLabel`s are `r paste(assessment$detID, collapse = ", ")`.  For simplicity, the term 'region' is used throughout to describe MIME regions.  OSPAR regions are always referred to as such.
</p>

<p class = "bot0">
The following sections:
</p>

<ul type = "disc">
<li>summarise the individual time series results</li>
<li>describe the methods used to combine the individual time series results to estimate regional trends and status</li>
<li>present the estimates of regional trends</li>
<li>present the estimates of regional status</li>
</ul>

<br>

### Summary of individual time series results {.tabset}

<p class = "startsmall">
A time series of `r detLabel` concentrations is assessed for status if:
</p>

<ul type = "disc">
<li>there is at least one year with data in the period `r periodStart` to `r periodEnd`</li>
<li>there are at least three years of data over the whole time series</li>
<li class = "last">a parametric model can be fitted to the data and used to estimate the mean concentration in the final monitoring year (or, occasionally, if a non-parametric test of status is applied)</li>
</ul>

<p class = "botsmall">
The conditions are more stringent for trends.  Specifically, a time series is assessed for trends if:
</p>

<ul type = "disc">
<li>there is at least one year with data in the period `r periodStart` to `r periodEnd`</li>
<li>there are at least five years of data over the whole time series</li>
<li class = "last">a parametric model can be fitted to the data and used to estimate the trend in mean concentrations</li>
</ul>

```{r include = FALSE, echo = FALSE}
nStatus <- with(assessment$timeSeries, sum(shape %in% c("filled_circle", "upward_triangle", "downward_triangle")))
nNonPar <- nStatus - with(assessment$timeSeries, sum(paramFit))
```

<p class = "botbig">
Note that all trend assessments for individual time series and most status assessments are based on the fit of a parametric model.  This is important because only the parametric results are passed into the regional assessments in the following sections.  Only `r nNonPar` of `r nStatus` `r detLabel` time series were assessed for status using a non-parametric test.
</p>

<p class = "botbig">
The first tab below shows all the monitoring stations where there are individual time series assessments of trend or status for `r detLabel`s. The magenta circles are stations where there is a trend assessment for at least one `r detLabel`, and the cyan cirles are stations where there are only status assessments.
</p>

<p class = "botsmall">
The next set of tabs show the trend and status assessments for each `r detLabel` in turn. The colours have the following meaning:
</p>

<ul class = "last">
```{asis, echo = isAC}
<li>blue: the mean concentration is significantly (p < 0.05) below the Background Assessment Concentration (BAC)</li>
<li>`r greenID`: the mean concentration is significantly (p < 0.05) below the `r detAC` (`r detAClabel`)</li>
<li>red: the mean concentration is not signficantly below the `r detAC` (`r detAClabel`)</li>
```

```{asis, echo = isBlack}
<li>black: there are no assessment criteria available</li>
```
</ul>

<p class = "botsmall">
And the shapes have the following meaning:
</p>

<ul type = "disc">
<li>downward triangle: the mean concentration is significantly decreasing (p < 0.05)</li>
<li>circle: there is no significant change in mean concentration (p > 0.05)</li>
<li class = "lastbig">upward triangle: the mean concentration is significantly increasing (p < 0.05)</li>
</ul>

<p class = "botsmall">
The final set of tabs give:
</p>

<ul class = "lastbig">
<li>the number of time series with upwards, downwards or no trend in each MiME and OSPAR region</li>
<li>the same information presented as proportions</li>
</ul>

<br>

#### Stations

```{r, echo = FALSE, message = FALSE, warning = FALSE, fig.width = 7, fig.height = 7}
plotStations(assessment, purpose = purpose, includeAll = TRUE)
```

<br>


```{r include = FALSE}
out <- NULL
for (detID in assessment$detID) {
  out <- c(out, knit_child('meta_analysis_bydet.Rmd'))
}
```


```{r, echo=FALSE, results="asis"}
cat(paste(out, collapse = '\n'))
```


#### Trend summary (numbers)

<p class = "startbig">
Number of time series with upwards, downwards or no trend by MIME and OSPAR region   
</p>

```{r echo = FALSE, message = FALSE, results = "asis"}
tab <- tableTimeSeries(assessment, "trend", asNumbers = TRUE, purpose = purpose)
tab <- xtable::xtable(tab)
print(tab, type = "html", include.rownames = FALSE) 
```

<br>

#### Trend summary (proportions)

<p class = "startbig">
Proportion of time series with upwards, downwards or no trend by MIME and OSPAR region  
</p>

```{r echo = FALSE, message = FALSE, results = "asis"}
tab <- tableTimeSeries(assessment, "trend", asNumbers = FALSE, purpose = purpose)
tab <- xtable::xtable(tab)
print(tab, type = "html", include.rownames = FALSE) 
```

<br>



### MIME regional assessment - methods 

<p class = "startbig">
Tabulating the number of time series with an upward or downward trend provides a quick summary of the individual time series results.  However, it does not provide an objective regional assessment of trends.  This section describes how the individual time series results can by synthesised in a meta-analysis to assess trends at the regional level.  However, in practice, there are only sufficient time series to do this for the Southern North Sea.  As no assessment criteria have yet been adopted for `r detLabel`s in `r params$compartment`, no formal regional assessment of status is possible.  However, instead, a meta-analysis can estimate the mean concentration profile of `r detLabel`s, again in the Southern North Sea.
</p>


#### Regional trend

<p class = "startsmall">
For a regional trend assessment, the trend in each time series is summarised by the estimated change in log concentration over the last twenty years (or shorter if the time series doesn't extend that far back). The mean trend for each organo-metal in the Southern North Sea is then estimated by fitting the following linear mixed model by restricted maximum likelihood:
</p>

<ul class = "last">
<li>response: trend (yearly change in log concentration)</li>
<li>fixed model: `r detLabel`</li>
<li>random model: station + trend estimation variation + residual variation</li>
</ul>

<p class = "botsmall">
The fixed model means that a separate regional trend is estimated for each `r detLabel`. The random model has three terms:
</p>

<ul type = "disc">
<li>station allows for variation in trend between stations common to all `r detLabel`s</li>
<li>trend estimation variation is the variance of the trend estimates from the individual time series analysis, assumed known and fixed</li>
<li class = "last">residual variation is the variation that cannot be explained by any of the fixed effects or the other random effects</li>
</ul>

<p class = "botlarge">
The meta-analysis is restricted to time series from monitoring stations that are classified as 'representative'.  'Baseline' stations (those with near pristine conditions or only subject to very diffuse inputs) are omitted because any trends there will likely be caused by different processes to those at representative stations.  Similarly 'impacted' stations (those close to a point source) are omitted because trends there will likely reflect changes due to the point source, and will be 'unrepresentative' at the regional level.  
</p>



#### Status

<p class = "startbig">
The profile of concentrations across `r detLabel`s can sometimes be useful for examining sources.  Profiles were examined by fitting the following linear mixed model by restricted maximum likelihood:

* response: mean log concentration 
* fixed model: `r detLabel`
* random model: station + estimation variation + residual variation

where, for each time series, the mean log concentration is the estimated mean log concentration in the final monitoring year and the estimation variance is the variance associated with this estimate.


#### Presentation

<p class = "startbig">
The results of both the trend and status regional assessments are back-transformed for presentation.  The regional trend is then interpreted as the percentage yearly change in concentration and the regional status is just the mean concentration. 
</p>

<br>


### Regional assessment - availability of time series {.tabset}

<p class = "startbig">
This section provides more detail on the number and geographic spread of the time series available for the meta-analyses.  The tabs show the stations with parametric trend and status assessments for each OSPAR region in turn, having excluded trend assessments at baseline and impacted stations.  The magenta circles are stations where there is a trend assessment for at least one `r detLabel`, and the cyan cirles are stations where there are only status assessments.
</p>

<p class = "botbig">
Note that, although there are no asssessment criteria available, the cyan circles show where there are time series with estimates of the mean concentration in the final montitoring year.  These can be used to construct mean concentration profiles.
</p>

<br>

```{r include = FALSE}
out <- NULL
for (regionID in as.character(1:4)) {
  out <- c(out, knit_child('meta_analysis_regional_plots.Rmd'))
}
```


```{r, echo=FALSE, results="asis"}
cat(paste(out, collapse = '\n'))
```




### Mean trends {.tabset}

<p class = "startsmall">
The first 2 tabs show:
</p>

<ul class = "last">
<li>Mean trend by `r detLabel`: mean trend estimates with pointwise `r ci` confidence limits by `r detLabel`</li>
<li>Individual time series: the `r plottype` estimates from the individual time series on which the meta-analysis is based.  The estimates from other regions (and excluded from the meta-analysis) are also shown.</li>
</ul>

<p class = "botsmall">
The symbols in these plots have the following interpretation:
</p>

<ul type = "disc">
<li>downward triangle: the mean concentration is significantly decreasing (p < 0.05)</li>
<li>circle: there is no change in mean concentration (p > 0.05)</li>
<li class = "lastbig">upward triangle: the mean concentration is significantly increasing (p < 0.05)</li>
</ul>

<p class = "startsmall">
The final tab shows the estimates of the mean trend by `r detLabel`.
</p>

<br>

#### Mean trends by `r detLabel`

```{r, echo = FALSE, message = FALSE, warning = FALSE, fig.width = 9, fig.height = 5}
plotTrendResults(assessment, "determinand")
```

<br>


#### Individual time series

```{r, echo = FALSE, message = FALSE, warning = FALSE, fig.width = 9, fig.height = 5}
plotTrend(assessment)
```

<br>

#### Estimates by `r detLabel`

<p class = "startsmall">
This table shows the estimates of the mean trend for each `r detLabel`:  
</p>

<ul type = "disc">
<li>trend: the estimated mean yearly change in log concentration across the region (multiplied by 100 for presentation)</li>
<li>se: the corresponding standard error</li>
<li>lower, upper: the corresponding pointwise 95% confidence limits</li>
<li>% yearly change: the estimated mean % yearly change in concentration across the region</li>
<li>%yc lower, %yc upper: the corresponding pointwise 95% confidence limits</li>
</ul>

<br>

```{r echo = FALSE, message = FALSE, results = "asis", eval = TRUE}
tab <- tableTrend(assessment, purpose = purpose, type = "determinand", detLabel = detLabel)
tab <- xtable::xtable(tab)
print(tab, type = "html", include.rownames = FALSE) 
```

<br>




### `r detLabel` profiles

The estimated `r detLabel` concentration profiles are shown below with pointwise 95% confidence limits.  The units are `r detUnit`.
```

```{r, echo = FALSE, message = FALSE, warning = FALSE, fig.width = 9, fig.height = 5}
plotStatusResults(assessment, "determinand", "CONC")
```


