Implementation of model routines from `mev` into `lax`

Question

Implementation of model routines from `mev` into `lax`

Closed this issue 5 years ago · 7 comments

Paul, what would be needed in the mev package for the classes mev_pp, mev_gpd, mev_rlarg and mev_gev to be supported by lax (outside of wrapper code to handle the arguments)?

These S3 classes have only two methods (print and plot), but I can extend this list.

Answer 1 · 2019-10-03T20:01:38.000Z

Sorry for being slow, Leo: I've just seen this today. I have pushed some code that supports mev_gev. The methods that you might want to add to mev are in mev_methods.R. Please feel free to copy these over to mev if you'd like. If you do this then I won't need to include them in lax. I'll do mev_gpd and mev_pp soon (it won't take long). I haven't yet incorporated any rlarg models, but this gives me a good incentive to do this.

(PS. I owe you an e-mail reply on another matter: I'll do this asap.)

Answer 2 · 2019-10-03T20:58:24.000Z

Thanks Paul (I had noted you were not following your own repo). I will add the methods for all four classes in the package for the next release and will send you the code in the meanwhile.

Answer 3 · 2019-10-03T22:50:01.000Z

Thank you. Yes, I need to follow my own repos! I'll sort that. If you could also make the following small change to mev::fit.pp() then it would enable me to deal add the class mev_pp. Please could you add to the returned list object the input data xdat? [I've just added mev_gpd. I'll look at rlarg asap.]

Answer 4 · 2019-10-04T00:37:04.000Z

I only keep threshold exceedances, so they are saved under exceedances rather than xdat, as for mev_gpd objects.

Answer 5 · 2019-10-04T11:20:12.000Z

Yes, that's fine for the GP model but if the argument cluster to alogLik() is non-NULL then for the PP model we need to know where in the original data the exceedances occurred. The non-exceedances contribute to the log-likelihood under the PP model and we need to be able to calculate the individual log-likelihood contributions from each of the observations, including the non-exceedances.

Answer 6 · 2019-10-04T18:57:39.000Z

So you only need xdat or you want something more? I could also return the indices of the exceedances rather than the original data if this is easier.

I should note that the mev_gpd method includes method="obre", which fits a robust estimator and the resulting covariance matrix is the sandwich variance estimator (from Hampel et al. (1986), eq. 4.2.13, p. 231).

I will also split the class methods for mev_gpdbayes to treat the approximate Bayesian inference methods differently, since they do not have vcov objects.

Answer 7 · 2019-10-07T07:27:37.000Z

Just xdat is fine: the code that I have use for other packages is set up to make the comparison with the threshold. I also need the argument npp (and perhaps np, although I haven't thought properly about np yet) to be returned. Sorry, I hadn't noticed that these are not returned at the moment.

Thanks for the heads up about method="obre". At the very least it I can use it as a check.

Cool, that makes sense. Many thanks.