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) fa.set = as.vector(unlist(FAset))

- 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) tombstone.info = predatorFAs[,1:4] predator.matrix = predatorFAs[,5:(ncol(predatorFAs))] npredators = nrow(predator.matrix)

- 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) prey.matrix = preyFAs[,-c(1,3)] # Selecting 5 prey species to include spec.red <-c("capelin", "herring", "mackerel", "pilchard", "sandlance") spec.red <- sort(spec.red) prey.red <- prey.matrix %>% 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.

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

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

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.

````r CalEst <- smufasa.est$Cal_Estimates

## Diet Estimates 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. ````r DietEst <- smufasa.est$Diet_Estimates

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

````r nll <- smufasa.est$nll

## Var_Epsilon This is the optimized diagonal values of the variance-covariance matrix of the errors. See reference in help file for details. ```r VarEps <- smufasa.est$Var_Epsilon

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.

Q = p.QFASA(predator.matrix=predator.matrix, prey.matrix=prey.matrix, 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:

DietEst.Q <- Q$'Diet Estimates'

*See* `p.QFASA`

*documnetation for further details.*

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

*See* `p.MUFASA`

*documnetation for further details.*

**Any scripts or data that you put into this service are public.**

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.