```
library(QFASA)
library(dplyr)
library(compositions)
```

Prior to starting make sure that:

- Fatty acid names in all files are the same (contain the exact same numbers/characters and punctuation)
- There are no fatty acids in the prey file that do not appear in the predator file and visa versa

- This is the list of FAs to be used in the modelling.
- The simplest alternative is to load a
`.csv`

file that contains a single column with a header row with the names of thee fatty acids listed below (see example file**“FAset.csv”**). - A more complicated alternative is to load a
`.csv`

file with the full set of FAs and then add code to subset the FAs you wish to use from that set ->*this alternative is useful if you are planning to test multiple sets*.

```
data(FAset)
= as.vector(unlist(FAset)) fa.set
```

- The FA signatures in the originating
`.csv`

file should be proportions summing to 1. - Each predator signature is a row with the FAs in columns (see example file
**“predatorFAs.csv”**). The FA signatures are subsetted for the chosen FA set (created above) and renormalized during the modelling so there is no need to subset and/or renormalize prior to loading the .csv file or running`p.SMUFASA`

BUT make sure that the same FAs appear in the predator and prey files. - Your predator FA
`.csv`

file can contain as much tombstone data columns as you like, you must extract the predator FA signatures as separate input in order to run in`p.SMUFASA`

. For example: in the code below, the predator`.csv`

file (“`predatorFAs.csv`

”) has 4 tombstone columns (SampleCode, AnimalCode, SampleGroup, Biopsy). Prior to running`p.SMUFASA`

, the tombstone (columns 1-4) and FA data (columns 5 onward) are each extracted from the original data frame. The FA data becomes the`predator.matrix`

(which is passed to`p.SMUFASA`

) and the tombstone data is retained so that it can be recombined with the model output later on.

```
data(predatorFAs)
= predatorFAs[,1:4]
tombstone.info = predatorFAs[,5:(ncol(predatorFAs))]
predator.matrix = nrow(predator.matrix) npredators
```

- The FA signatures in the
`.csv`

file should be proportions summing to 1. - The prey file should contain all of the individual FA signatures of the prey and their lipid contents (where appropriate).
- If you want to only include a subset of prey species, you must extract it prior to input (see code below).

```
data(preyFAs)
= preyFAs[,-c(1,3)]
prey.matrix
# Selecting 5 prey species to include
<-c("capelin", "herring", "mackerel", "pilchard", "sandlance")
spec.red <- sort(spec.red)
spec.red <- prey.matrix %>%
prey.red filter(Species %in% spec.red)
```

- Mean lipid content by species group is calculated from the full prey file using the species group as a summary variable (see code below).
**Note**: If no lipid content correction is going to be applied, then a vector of 1s of length equal to the number of species groups is used as the vector instead. I.e.`FC - rep(1,nrow(prey.matrix))`

.- If you’ve decided on a subset of species, you must extract them from the mean lipid content vector as well.

```
= preyFAs[,c(2,3)]
FC = FC %>%
FC arrange(Species)
= tapply(FC$lipid,FC$Species,mean,na.rm=TRUE)
FC.vec <- FC.vec[spec.red] FC.red
```

`<- p.SMUFASA(predator.matrix, prey.red, FC.red, fa.set) smufasa.est `

The MUFASA output is a list with 4 components:

- Cal_Estimates
- Diet_Estimates
- nll
- Var_Epsilon

A vector of length equal to the number of FAs used and whose sum is the total number of FAs used. Thos is a set of calibration coefficients common to all predators used.

`<- smufasa.est$Cal_Estimates CalEst `

The diet estimate vector returned by `p.SMUFASA`

represents an overall common diet for all predators in `predator.matrix`

. **Note**: If you wish to obtain diet estimates for each individual predator see the steps below.

`<- smufasa.est$Diet_Estimates DietEst `

This is a vector of the negative log likelihood values at each iteration of the optimizer.

`<- smufasa.est$nll nll `

This is the optimized diagonal values of the variance-covariance matrix of the errors. See reference in help file for details.

`<- smufasa.est$Var_Epsilon VarEps `

Once a vector of calibration coefficients is obtained via `p.SMUFASA`

you can pass this vector to `p.QFASA`

or `p.MUFASA`

to obtain individual diet estimates.

```
= p.QFASA(predator.matrix=predator.matrix, prey.matrix=prey.matrix,
Q cal.mat=CalEst, dist.meas=2, gamma=1, FC=FC.red,
start.val=rep(1,nrow(prey.red)), ext.fa=fa.set)
```

QFASA Diet Estimates:

`<- Q$'Diet Estimates' DietEst.Q `

*See* `p.QFASA`

*documnetation for further details.*

```
<- p.SMUFASA(predator.matrix, prey.red, cal.est, FC.red, fa.set)
M <- M$Diet_Estimates DietEst.M
```

*See* `p.MUFASA`

*documnetation for further details.*