[R] Peak Over Threshold values

Clint Bowman clint at ecy.wa.gov
Wed May 26 20:52:51 CEST 2010


How about?

hi.rle<-rle(walevel>5.79)
lo.rle<-rle(walevel<5.36)
plot(walevel)
abline(h=5.8,col=2,lty=3)
abline(h=5.35,col=3,lty=3)
hi.lo.rle<-sort(c(cumsum(hi.rle$lengths),cumsum(lo.rle$lengths)))
abline(v=hi.lo.rle)

You can use the $values from the rle to sort things out.  Probably 
want to ignore the end point at length(walevel).

-- 
Clint Bowman			INTERNET:	clint at ecy.wa.gov
Air Quality Modeler		INTERNET:	clint at math.utah.edu
Department of Ecology		VOICE:		(360) 407-6815
PO Box 47600			FAX:		(360) 407-7534
Olympia, WA 98504-7600

On Wed, 26 May 2010, Vito Muggeo (UniPa) wrote:

> dear Tonja,
> By plotting your data
>
> plot(df)
>
> it seems to me that you are looking for a piecewise linear relationships. If 
> this is the case, have a look to the package segmented. You have to specify 
> or not the number and the starting values for the breakpoints
>
> library(segmented)
> olm<-lm(walevel~day)
>
> #specify number and starting values for the breakpoints..
> oseg<-segmented(olm, seg.Z=~day, psi=c(20,50,80))
> plot(oseg,add=TRUE,col=2)
> oseg$psi
>
> #or you can avoid to specify starting values for psi
> oseg<-segmented(olm, seg.Z=~day, 
> psi=NA,control=seg.control(stop.if.error=FALSE,K=30))
>
> plot(oseg,add=TRUE,col=2)
> oseg$psi
>
>
> best,
> vito
>
>
> Tonja Krueger ha scritto:
>>     Dear List
>>
>>     I hope you can help me: I’ve got a dataframe (df) within which I am
>>     looking
>>     for Peak Over Threshold values as well as the length of the events. An
>>     event
>>     starts when walevel equals 5.8 and it should end when walevel equals
>>     the
>>     lower threshold value (5.35).
>>
>>     I tried “clusters (…)” from “evd package”, and varied r (see example)
>>     but it
>>     did not work for all events (again see example).
>>
>>     walevel <- c(5.75, 5.75, 5.75, 5.76, 5.80, 5.82, 5.85, 5.86, 5.91,
>>     5.91,
>>     5.88, 5.92, 5.99, 6.02, 6.00, 6.08, 6.11, 6.10, 6.10, 6.11, 6.14, 6.12,
>>     6.15, 6.17, 6.15, 6.08, 6.01, 5.95, 5.89, 5.83, 5.77, 5.72, 5.70, 5.65,
>>     5.59, 5.51, 5.43, 5.33, 5.30, 5.25, 5.22, 5.21, 5.19, 5.19, 5.19, 5.13,
>>     5.15, 5.07, 5.04, 5.02, 4.99, 5.05, 5.07, 5.12, 5.17, 5.22, 5.32, 5.29,
>>     5.33, 5.33, 5.36, 5.37, 5.42, 5.43, 5.48, 5.55, 5.57, 5.66, 5.68, 5.72,
>>     5.77, 5.79, 5.81, 5.80, 5.80, 5.77, 5.72, 5.70, 5.66, 5.62, 5.62, 5.61,
>>     5.59, 5.57, 5.51, 5.47, 5.41, 5.39, 5.40, 5.40, 5.42, 5.45, 5.47, 5.50,
>>     5.50, 5.49, 5.43, 5.39, 5.33, 5.26)
>>
>>     day <- c(1:100)
>>
>>     df <- data.frame(day,walevel)
>>
>>     library(evd)
>>     clusters(df$walevel, u = 5.80, r = 1, ulow = 5.35, cmax = T, plot = T)
>>     clusters(df$walevel, u = 5.80, r = 50, ulow = 5.35, cmax = T, plot = T)
>>
>>     What have I done wrong?
>>
>>     Tonja
>>  ______________________________________________
>>  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.
>> 
>
>


More information about the R-help mailing list