Here is a “solution” for practical 3. As with any data analysis, there is no correct answer, but this shows how I would approach this analysis. The analysis here is conditional on selecting a detection function in the previous exercises; I’ve chosen `df_hn`

and `df_hr_ss_size`

that we saved previously (see detection function solutions).

Much of the text below is as in the exercise itself, so it should be relatively easy to navigate.

Additional text and code is highlighted using boxes like this.

`library(Distance)`

```
## Loading required package: mrds
## This is mrds 2.1.14
## Built: R 3.2.0; ; 2015-07-30 10:07:19 UTC; unix
```

`library(dsm)`

```
## Loading required package: mgcv
## Loading required package: nlme
## This is mgcv 1.8-7. For overview type 'help("mgcv-package")'.
## This is dsm 2.2.11
## Built: R 3.2.2; ; 2015-10-23 20:20:41 UTC; unix
```

```
library(ggplot2)
library(knitr)
```

Loading the `RData`

files where we saved our results:

```
load("sperm-data.RData")
load("df-models.RData")
```

Before we fit a model using `dsm()`

we must first remove the observations from the spatial data that we excluded when we fitted the detection function – those observations at distances greater than the truncation.

`obs <- obs[obs$distance <= df_hn$ddf$meta.data$width,]`

Here we’ve used the value of the truncation stored in the detection function object, but we could also use the numeric value (which we can also find by checking the model’s `summary()`

).

Also note that if you want to fit DSMs using detection functions with different truncation distances, then you’ll need to reload the `sperm-data.RData`

and do the truncation again for that detection function.

Using the data that we’ve saved so far, we can build a call to the `dsm()`

function and fit out first density surface model. Here we’re only going to look at models that include spatial smooths.

Let’s start with a very simple model – a bivariate smooth of `x`

and `y`

:

```
dsm_nb_xy <- dsm(count~s(x,y),
ddf.obj=df_hn, segment.data = segs, observation.data=obs,
family=nb(), method="REML")
```

Note again that we try to have informative model object names so that we can work out what the main features of the model were from its name alone.

We can look at a `summary()`

of this model. Look through the summary output and try to pick out the important information based on what we’ve talked about in the lectures so far.

`summary(dsm_nb_xy)`

```
##
## Family: Negative Binomial(0.105)
## Link function: log
##
## Formula:
## count ~ s(x, y) + offset(off.set)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -20.7009 0.2538 -81.56 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(x,y) 17.95 22.23 75.89 6.27e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.0879 Deviance explained = 40.5%
## -REML = 392.65 Scale est. = 1 n = 949
```

Important things to look at here include the

`EDF`

column, the `R-sq.(adj)`

and `Deviance explained`

. Note that the significance of the terms isn’t that useful at the moment, since we only have one smooth in the model (though if this wasn’t significant we might be worried about including that variable!).
As discussed in the lectures, the `plot`

output is not terribly useful for bivariate smooths like these. We’ll use `vis.gam()`

to visualise the smooth instead:

`vis.gam(dsm_nb_xy, view=c("x","y"), plot.type="contour", too.far=0.1, main="s(x,y) (link scale)", asp=1)`

Notes:

- The plot is on the scale of the link function, the offset is not taken into account – the contour values do not represent abundance, just the “influence” of the smooth.
- We set
`view=c("x","y")`

to display the smooths for`x`

and`y`

(we can choose any two variables in our model to display like this) `plot.type="contour"`

gives this “flat” plot, set`plot.type="persp"`

for a “perspective” plot, in 3D.- The
`too.far=0.1`

argument displays the values of the smooth not “too far” from the data (try changing this value to see what happens. `asp=1`

ensures that the aspect ratio of the plot is 1, making the pixels square.- Read the
`?vis.gam`

manual page for more information on the plotting options.

We can use the `gam.check()`

and `rqgam.check`

functions to check the model.

`gam.check(dsm_nb_xy)`

```
##
## Method: REML Optimizer: outer newton
## full convergence after 7 iterations.
## Gradient range [-1.910812e-06,2.784947e-06]
## (score 392.646 & scale 1).
## Hessian positive definite, eigenvalue range [2.157928,29.21001].
## Model rank = 30 / 30
##
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
##
## k' edf k-index p-value
## s(x,y) 29.000 17.953 0.532 0
```

`rqgam.check(dsm_nb_xy)`