Peter Maclean
2011-Jun-30 04:16 UTC
[R] Saving fExtremes estimates and k-block return level with confidence intervals.
I am estimating a large model by groups. How do you save the results and?returns the associated quantiles? For this example I need a data frame n?? ?xi??????? mu????????beta 1?? 0.1033614? 2.5389580 0.9092611 2? ?0.3401922? 0.5192882 1.5290615 3?? 0.5130798? 0.5668308 1.2105666 I also want to apply gevrlevelPlot() for each "n" or group. ? #Example n <- c(1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3,3) y <- c(2,3,2,3,4,5,6,1,0,0,0,6, 2, 1, 0, 0,9,3) z <- as.data.frame(cbind(n,y)) colnames(z) <- c("n","y") library(fExtremes) z <- split(z, z$n) res2 <-lapply(z, function(x){ ?????????????? m <- as.numeric(x$y) ?????????????? gevFit(m, block = 1, type = c("pwm")) ??????????????? })> res2$`1` Title: ?GEV Parameter Estimation Call: ?gevFit(x = m, block = 1, type = c("pwm")) Estimation Type: ? gev pwm Estimated Parameters: ?????? xi??????? mu????? beta 0.1033614 2.5389580 0.9092611 Description ? Wed Jun 29 23:07:48 2011 $`2` Title: ?GEV Parameter Estimation Call: ?gevFit(x = m, block = 1, type = c("pwm")) Estimation Type: ? gev pwm Estimated Parameters: ?????? xi??????? mu????? beta 0.3401922 0.5192882 1.5290615 Description ? Wed Jun 29 23:07:48 2011 $`3` Title: ?GEV Parameter Estimation Call: ?gevFit(x = m, block = 1, type = c("pwm")) Estimation Type: ? gev pwm Estimated Parameters: ?????? xi??????? mu????? beta 0.5130798 0.5668308 1.2105666 Description ? Wed Jun 29 23:07:48 2011?
Jorge Ivan Velez
2011-Jun-30 07:37 UTC
[R] Saving fExtremes estimates and k-block return level with confidence intervals.
Hi Peter, Try data.frame(n = names(res2), t(sapply(res2, function(l) l@fit$par.ests))) for the first part. HTH, Jorge On Thu, Jun 30, 2011 at 12:16 AM, Peter Maclean <> wrote:> I am estimating a large model by groups. How do you save the results > and returns > the associated quantiles? > For this example I need a data frame > n xi mu beta > 1 0.1033614 2.5389580 0.9092611 > 2 0.3401922 0.5192882 1.5290615 > 3 0.5130798 0.5668308 1.2105666 > I also want to apply gevrlevelPlot() for each "n" or group. > > #Example > n <- c(1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3,3) > y <- c(2,3,2,3,4,5,6,1,0,0,0,6, 2, 1, 0, 0,9,3) > z <- as.data.frame(cbind(n,y)) > colnames(z) <- c("n","y") > library(fExtremes) > z <- split(z, z$n) > res2 <-lapply(z, function(x){ > m <- as.numeric(x$y) > gevFit(m, block = 1, type = c("pwm")) > }) > > res2 > $`1` > Title: > GEV Parameter Estimation > Call: > gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > gev pwm > Estimated Parameters: > xi mu beta > 0.1033614 2.5389580 0.9092611 > Description > Wed Jun 29 23:07:48 2011 > > $`2` > Title: > GEV Parameter Estimation > Call: > gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > gev pwm > Estimated Parameters: > xi mu beta > 0.3401922 0.5192882 1.5290615 > Description > Wed Jun 29 23:07:48 2011 > > $`3` > Title: > GEV Parameter Estimation > Call: > gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > gev pwm > Estimated Parameters: > xi mu beta > 0.5130798 0.5668308 1.2105666 > Description > Wed Jun 29 23:07:48 2011 > > > ______________________________________________ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. >[[alternative HTML version deleted]]
Peter Maclean
2011-Jul-06 08:29 UTC
[R] Saving fExtremes estimates and k-block return level with confidence intervals.
Hi: I am trying to compare the results of "evir" and "fExtreme" packages. I?could not figure out how to save the "evir" package results.?Also, how to pass the results to?"fExtreme" function "gevrlevelPlot" and?"evir" function "rlevel.gev" to get the return levels. I just need the?values and not graphs. ? #Example library(fExtremes) library(evir) require(plyr) y<- data.frame(rgev(300, xi = 0.1, mu = .5, sigma = 1.6)) colnames(y) <- c("y") n <- data.frame(n = rep(c(1,2,3), each = 100)) z <- cbind(n,y) y <- z$y # Model for grouped data ##fExtremes package z1 <- split(z,z$n) rgf <- lapply(z1, function(x){ ????????????? m <- as.numeric(x$y) ????????????? gevFit(m, block = 2, type = "mle") ?}) #Save results resf<- ldply(rgf, function(x) x at fit$par.ests) #Qs: How to transfer rge object to "gevrlevelPlot" function to get the values? ? ##evir package rge<- by(z,z[,"n"], function(x) gev(x$y, 2, method = "BFGS", control =list(maxit = 10000)))? ? #Qs:How to save par.ests and?par.ses? #Qs:How to transfer rge object to "rlevel.gev" function to get the values??????????? ? #Model for single vector rlf <- gevFit(y, block = 2, type = "mle") rlfp <- gevrlevelPlot(rlf, kBlocks = 2) rlfp rle <- gev(y, 2, method = "BFGS", control =list(maxit = 10000)) rlep<- rlevel.gev(rle, k.blocks = 2, add = FALSE) rlep ----- Original Message ---- From: Dennis Murphy <djmuser at gmail.com> To: Peter Maclean <pmaclean2011 at yahoo.com> Sent: Thu, June 30, 2011 3:03:05 PM Subject: Re: [R] Saving fExtremes estimates and k-block return level with confidence intervals. Hi: The plyr package can help you out as well. As far as the estimates go, library(plyr)> ldply(res2, function(x) x at fit$par.ests)? .id? ? ? ? xi? ? ? ? mu? ? ? beta 1? 1 0.1033614 2.5389580 0.9092611 2? 2 0.3401922 0.5192882 1.5290615 3? 3 0.5130798 0.5668308 1.2105666 You could also look into the l_ply() function, which takes a list as input and outputs nothing; however, it is usually called if one want to make a series of similar plots. You need to write a function that takes a generic list component and outputs a plot. Here's a fairly trivial example: testdf <- data.frame(gp = rep(LETTERS[1:3], each = 10), x = 1:10, ? ? ? ? ? ? ? ? ? ? ? y = rnorm(30)) l_ply(split(testdf, testdf$gp), function(df) plot(y ~ x, data = df)) Alternatively, # Observe that a data frame is used as input pfun <- function(df) plot(y ~ x, data = df) l_ply(split(testdf, testdf$gp), pfun) In this case, l_ply() plots y vs. x in each of the three subgroups of testdf separately. It appears you want to do a similar thing, but with the list of model outputs instead; in your case, the function you write should take a list as its sole argument. Any slots or components that are accessed are then relative to the input list object - see the anonymous function I wrote for the output data frame above as an illustration. Since it appears that gevrlevelPlot() also outputs a vector, you may want to write a function that returns the output vector as a data frame and call ldply() instead. I don't know enough about the fExtremes package to write this myself, but it should be possible to write a function for a generic model that generates both a plot and an output data frame. The plyr package is pretty good at this sort of thing. HTH, Dennis On Wed, Jun 29, 2011 at 9:16 PM, Peter Maclean <pmaclean2011 at yahoo.com> wrote:> I am estimating a large model by groups. How do you save the results >and?returns > the associated quantiles? > For this example I need a data frame > n?? ?xi??????? mu????????beta > 1?? 0.1033614? 2.5389580 0.9092611 > 2? ?0.3401922? 0.5192882 1.5290615 > 3?? 0.5130798? 0.5668308 1.2105666 > I also want to apply gevrlevelPlot() for each "n" or group. > > #Example > n <- c(1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3,3) > y <- c(2,3,2,3,4,5,6,1,0,0,0,6, 2, 1, 0, 0,9,3) > z <- as.data.frame(cbind(n,y)) > colnames(z) <- c("n","y") > library(fExtremes) > z <- split(z, z$n) > res2 <-lapply(z, function(x){ > ?????????????? m <- as.numeric(x$y) > ?????????????? gevFit(m, block = 1, type = c("pwm")) > ??????????????? }) >> res2 > $`1` > Title: > ?GEV Parameter Estimation > Call: > ?gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > ? gev pwm > Estimated Parameters: > ?????? xi??????? mu????? beta > 0.1033614 2.5389580 0.9092611 > Description > ? Wed Jun 29 23:07:48 2011 > > $`2` > Title: > ?GEV Parameter Estimation > Call: > ?gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > ? gev pwm > Estimated Parameters: > ?????? xi??????? mu????? beta > 0.3401922 0.5192882 1.5290615 > Description > ? Wed Jun 29 23:07:48 2011 > > $`3` > Title: > ?GEV Parameter Estimation > Call: > ?gevFit(x = m, block = 1, type = c("pwm")) > Estimation Type: > ? gev pwm > Estimated Parameters: > ?????? xi??????? mu????? beta > 0.5130798 0.5668308 1.2105666 > Description > ? Wed Jun 29 23:07:48 2011 > > > ______________________________________________ > R-help at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. >