[FLR-list] fitting a custom SRModel with a covariate

Iago MOSQUEIRA (JRC) iago.mosqueira at jrc.ec.europa.eu
Thu Apr 7 16:46:08 CEST 2016


Hi,

trying the example I found a bug in FLCore. This just shows that no much 
use has been made of this model ...

Can you install from github (i.e. do you have Rtools installed if in 
windows?)? If so do

library(devtools)

install_github('flr/FLCore')

and then

library(FLCore)

data("ple4")

ple4SR <- as.FLSR(ple4)

env1 <- c(-0.81, -1.73, -0.09, 0.42, -1.01, -1.06, -0.79, -0.42, -0.74, 
-0.3, 0.25, 0.09, -0.18, -0.29, 0.13, 0.09, 0.48, 0.33, 0.41, 0.4, 
-0.21, -0.64, -0.97, -0.02, -0.19, 0.28, -0.03, -0.19, -0.47, -0.97, 
-0.88, 0.23, 0.82, 0.92, 0.17, 0.67, -0.07, 0.32, 0.9, -0.64, 0.97, 
0.54, 1.18, 0.67, 0.44, 1.29, 1.54, 1.22, 0.47, 1.04, 1.31, 1.06)

covar(ple4SR) <- FLQuants(covar=FLQuant(env1, 
dimnames=dimnames(ssb(ple4SR))))

model(ple4SR) <- rickerCa

ple4SR <- fmle(ple4SR)

plot(ple4SR)

This now fits, but I am not so sure about the values. Does this covar 
make sense for these data? Or is it just for testing?

Let me know how it goes. It you cannot install from source, I will 
compile a new package tomorrow.

Cheers,


Iago



On 07/04/16 17:03, Marc Taylor wrote:
> Thanks Laurie. I would greatly appreciate the example if you are able to
> find it.
>
> Cheers,
> Marc
>
>
> Achtung: Das Thünen-Institut hat die Domain gewechselt. Bitte ändern Sie
> meine Mailadresse in Ihrem Adressbuch!
> Notice: Thünen Institute has changed its domain. Please change my email
> address in your address book!
>
> Dr. Marc Taylor
> Marine Lebende Resourcen / Marine Living Resources
> Thünen-Institut für Seefischerei / Thünen Institute of Sea Fisheries
> Palmaille 9
> 22767 Hamburg, Germany
>
> Tel:  +49 40 38905-143
> Mail: marc.taylor at thuenen.de
> Web:  www.ti.bund.de
>
> On 4/7/2016 9:49 AM, laurie wrote:
>> see, which is in FLCore
>>
>>
>> rickerCa <- function() {
>>   logl <- function(a, b, c, rec, ssb, covar)
>>     loglAR1(log(rec), log(a * (1 - c * covar) * ssb * exp(-b * ssb)))
>>
>>   initial <- structure(function(rec, ssb) {
>>         # The function to provide initial values
>>     res  <-coefficients(lm(c(log(rec/ssb))~c(ssb)))
>>     return(FLPar(a=max(exp(res[1])), b=-max(res[2]), c=1))},
>>
>>   # lower and upper limits for optim()
>>     lower=rep(-Inf, 3),
>>     upper=rep( Inf, 3))
>>
>>     model  <- rec ~ a * (1 - c * covar) * ssb * exp(-b * ssb)
>>     return(list(logl=logl, model=model, initial=initial))
>> }
>>
>>
>> There is an example using it somewhere, I will see if I can find it.
>>
>> Laurie
>>
>> On 07/04/16 09:37, Marc Taylor wrote:
>>> Dear FLR-Listers,
>>>
>>> I was hoping someone might help me figure out a way to create a
>>> custom SRModel that uses a covariate in addition to ssb and rec.
>>>
>>> From the documentation of SRModels (i.e. ?SRModels), reference is
>>> made to several models that would be of interest to me (e.g.
>>> 'ricker.c.a', 'ricker.c.b', etc.), but these do not appear to be
>>> presently included in FLCore.  I made an attempt (below) to define a
>>> custom model, but am getting errors when fitting with fmle().
>>>
>>> If anyone has any experience with such models, I would greatly
>>> appreciate any advise.
>>>
>>> Cheers,
>>> Marc
>>>
>>> *Example script*:
>>>
>>> # packages
>>> ----------------------------------------------------------------
>>> library(FLCore) # FLCore_2.5.20160107
>>>
>>>
>>> # load data
>>> ---------------------------------------------------------------
>>>
>>> data("ple4")
>>>
>>>
>>> # make srr obbject and add covariate
>>> --------------------------------------
>>>
>>> ple4Sr <- as.FLSR(ple4)
>>>
>>> env1 <- c(-0.81, -1.73, -0.09, 0.42, -1.01, -1.06, -0.79, -0.42, -0.74,
>>> -0.3, 0.25, 0.09, -0.18, -0.29, 0.13, 0.09, 0.48, 0.33, 0.41,
>>> 0.4, -0.21, -0.64, -0.97, -0.02, -0.19, 0.28, -0.03, -0.19, -0.47,
>>> -0.97, -0.88, 0.23, 0.82, 0.92, 0.17, 0.67, -0.07, 0.32, 0.9,
>>> -0.64, 0.97, 0.54, 1.18, 0.67, 0.44, 1.29, 1.54, 1.22, 0.47,
>>> 1.04, 1.31, 1.06)
>>>
>>> ple4Sr at covar <- FLQuants(list(env1=FLQuant(env1, dim=dim(ple4Sr at ssb),
>>> dimnames=dimnames(ple4Sr at ssb))))
>>> ple4Sr at covar$env1
>>>
>>>
>>> # attempt with existing SRModel fails
>>> model(ple4Sr) <- "ricker.c.a"
>>> # Error in do.call(value, list()) : could not find function "ricker.c.a"
>>>
>>>
>>>
>>> # make new srr function
>>> ---------------------------------------------------
>>>
>>> # adapt ricker
>>> ricker <- function(){
>>>     logl <- function(a, b, rec, ssb) loglAR1(log(rec), log(a *
>>>         ssb * exp(-b * ssb)))
>>>     initial <- structure(function(rec, ssb) {
>>>         res <- coefficients(lm(log(c(rec)/c(ssb)) ~ c(ssb)))
>>>         return(FLPar(a = max(exp(res[1])), b = -max(res[2])))
>>>     }, lower = rep(-Inf, 2), upper = rep(Inf, 2))
>>>     model <- rec ~ a * ssb * exp(-b * ssb)
>>>     return(list(logl = logl, model = model, initial = initial))
>>> }
>>>
>>> # to this (incl. env. control term - 'env1')
>>> ri.ec <- function(){
>>>     logl <- function(a, b, e1, rec, ssb, env1) loglAR1(log(rec), log(a *
>>>         ssb * exp(-b * ssb) * exp(e1 * env1)))
>>>     initial <- structure(function(rec, ssb, env1) {
>>>         res <- coefficients(lm(log(c(rec)/c(ssb)) ~ c(ssb) + c(env1)))
>>>         return(FLPar(a = max(exp(res[1])), b = -max(res[2], env1 =
>>> max(res[3]))))
>>>     }, lower = rep(-Inf, 3), upper = rep(Inf, 3))
>>>     model <- rec ~ a * ssb * exp(-b * ssb) * exp(e1 * env1)
>>>     return(list(logl = logl, model = model, initial = initial))
>>> }
>>>
>>>
>>> # fit SRR
>>> -----------------------------------------------------------------
>>> model(ple4Sr) <- "ri.ec"
>>> ple4Sr <- fmle(ple4Sr)
>>> # Error in loglAR1(log(rec), log(a * ssb * exp(-b * ssb) * exp(e1 *
>>> env1))) :
>>> #   error in evaluating the argument 'hat' in selecting a method for
>>> function 'loglAR1': Error: argument "e1" is missing, with no default
>>> summary(ple4Sr)
>>> plot(ple4Sr)
>>> --
>>> Achtung: Das Thünen-Institut hat die Domain gewechselt. Bitte ändern
>>> Sie meine Mailadresse in Ihrem Adressbuch!
>>> Notice: Thünen Institute has changed its domain. Please change my
>>> email address in your address book!
>>>
>>> Dr. Marc Taylor
>>> Marine Lebende Resourcen / Marine Living Resources
>>> Thünen-Institut für Seefischerei / Thünen Institute of Sea Fisheries
>>> Palmaille 9
>>> 22767 Hamburg, Germany
>>>
>>> Tel:  +49 40 38905-143
>>> Mail:marc.taylor at thuenen.de
>>> Web:www.ti.bund.de
>>>
>>>
>>> _______________________________________________
>>> flr-list mailing list
>>> flr-list at flr-project.org
>>> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/flr-list
>>
>>
>>
>> _______________________________________________
>> flr-list mailing list
>> flr-list at flr-project.org
>> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/flr-list
>
>
>
>
> _______________________________________________
> flr-list mailing list
> flr-list at flr-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/flr-list
>


-- 
Dr Iago Mosqueira

European Commission – Joint Research Center
Institute for the Protection and Security of the Citizen (IPSC)
Maritime Affairs Unit G03
TP 051, Via Enrico Fermi 2749
I-21027 Ispra (VA), Italy

Office : +39 0332 785413
Fax: +39 0332 789658

iago.mosqueira at jrc.ec.europa.eu

http://fishreg.jrc.ec.europa.eu/home


More information about the flr-list mailing list