
```{r setup44, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning=FALSE,message = FALSE)
```

# Arne data analysis part one


```{r}
library(mapview)
library(leaflet.extras)
library(tidyverse)
library(mgcv)
library(aqm)
library(sf)
```


In 2017, 2018 and 2019 students on this unit conducted field work at Arne in order to collect some original data for analysis.

From the perspective of the ILOs of the unit the objective of the study was to put into practice some of the techniques for analysing data using both spatially explicit and no spatial techniques.

The questions being addressed by the exercise revolved around the management of lowland heathland, but they also have some broader ecological relevance.

Scots pine (*Pinus sylvestris*) is native to the UK. There is evidence that the species was once quite widely distributed Southern England and may have been a common species in the neolithic landscape [@whitehouse2010] Habitats in the UK have been influenced by human activity throughout the holocene, in a wide range of different ways. Scots pine appears to have been extirpated from most England by human over exploitation and clearance for agriculture by around 6500 years before present leading to the widespread development of open heathland by around 3000 years before present[@groves2012]. The precise cause and date of the loss is not entirely clear. Pines were reintroduced for timber in the 17th century and those found at Arne probably have an originally central European genetic provenance. They are therefore often regarded as an introduced species. The pine forests in Dorset are all the result of deliberate management for timber. Lowland heath is a priority habitat for conservation and management at Arne is focussed on the restoration of open heathland, with some scattered trees. A challenge for management arises as a result of natural regeneration of pines in the open heath and difficulties in retaining previously forested areas as heathland after felling.

In order to keep the heathland open the RSPB use periodic "pine pulling" campaigns. There are even days on which you can visit the reserve to "pull a pine for Christmas". Part two of the exercise will look at the pine regeneration on the heath in more detail.

Pines establishing on heathland may show less vigorous growth than those establishing after felling. A possible explanation for this may be the presence of mycorrhizal fungi that mediate invasion [@collier2009]. 

During the 2019 field trip the students looked at the differences between the growth of young pines in two areas at Arne. One area was heathland with scattered young pines that regenerated from seed sources on the heath. The second site was once a pine forest which was cleared in 2012 in order to restore heathland. Many pines have regrown since from the seedbank in the soil. 

The aim was to quantify any detectable difference using a study that could be designed to be conducted in less than two hours. As it was raining heavily the study was restricted to two small areas. It is worth pointing out that this causes some difficulties for formal inference due to spatial autocorrelation. You can ignore this for the purposes of the exercise.

The diameter of young regenerating pines 3cm above the ground was measured together with their height. An estimate of age was also made by counting the number of whorls. GPS way points were recorded and matched to the data. Pines that are growing vigorously would be expected to be taller for the same diameter as pines that are growing in less favourable conditions. This is typical metric used by foresters to evaluate site quality and should apply to the small pines.

DiMM : Diameter in mm  
HeCM : Height in cm



Run the following code to load the data. You will need the sf and mapview packages loaded.

```{r}
data(arne_pine_hts)
pine_hts<-st_transform(pine_hts,27700) # British national grid
pine_hts$Site<-as.factor(pine_hts$Site) # Convert the numerical value of the site to a factor
 pine_hts %>% st_drop_geometry() -> d ## Take out geometry and make a simple data frame called d
 
```

You can now look at the positions of the data points using mapview.

```{r}
mapview(pine_hts,z="Site") ->map
map@map%>% addMiniMap(position = "bottomleft",zoomLevelOffset = -4, toggleDisplay = TRUE, minimized = TRUE)
```

Pines growing at site 2 had regenerated following clearance. The clearance of the pines even shows up on a global classification of Landsat derived satellite imagery [@hansen2013]

```{r}
mapview(pine_hts,z="Site") %>% giscourse::hansen_wms()
```



## Regression

One way of looking at differences in growth between the sites may be to first find a relationship between diamater and height and then obtain the residuals (there are other options).

You can try the line and spline approach.

```{r, eval=FALSE}
d %>% ggplot(aes(x=DiMM, y=HeCM)) ->g0
g0 + geom_point() + geom_smooth(method="gam", formula =y~s(x)) + 
  geom_smooth(method="lm", colour="red")
```

Is a linear model a suitable approximation over this data range?

We can take the residuals of the models and then apply the techniques you have learned before.

```{r}
modgam<-gam(data=d,HeCM~s(DiMM)) # The gam model
modlm<-lm(data=d,HeCM~DiMM)   # The linear model

d$residuals<-residuals(modlm)  # Take the residuals of the preferred mode: Swap for modgam if you want to. 
```

Now think about how you can look at the differences between sites.

## Analysis of covariance

An alternative to this would be to look at the differences between the slopes of a linear model fitted to each site.


```{r, eval=FALSE}
d %>% ggplot(aes(x=DiMM, y=HeCM, colour=Site)) ->g0
g0 +  geom_point() +
  geom_smooth(method="lm") 
```

Alternatively you could use a facet wrap, which is better if the figure becomes cluttered.

```{r, eval=FALSE}
d %>% ggplot(aes(x=DiMM, y=HeCM)) ->g0
g0 +  geom_point() +
  geom_smooth(method="lm", colour="red") + facet_wrap(~Site)
```


They do look different. You can check this using a technique called **analysis of covariance**. If the interaction term is significant then there is a detectable difference. The theory behind this is rather advanced, but it is worth knowing that this technique is available.  To run an analysis of covariance you combine a numerical variable with a factor.

```{r, eval=FALSE}
mod<-lm(data=d,HeCM~DiMM*Site) 
anova(mod)
```


Notice that the interaction term DiMM:Site is highly significant. This implies that the two slopes do differ. This is observable in the figures, so the analysis formalises this. 

Now we can find a confidence interval for the differences between the two slopes.

```{r, eval=FALSE}
confint(mod)
```
The CI for DiMM shown here represents the confidence interval for the slope not mentioned (i.e. Site 1). The DiMM:Site2 sonfidence interval represents the CI for the difference between the slopes for the line fitted for Site1 and Site2.


## Exercise

Complete an analysis of these data using the code shown here and additional code adapted from previous exercises. As this is evaluated no direct guidance regarding how to apply the methods will be provided, as this is for you to decide.


