[R] help with script to get starting date of blooms

Jim Lemon jim at bitwrit.com.au
Sat Jul 19 13:25:01 CEST 2014


On Fri, 18 Jul 2014 10:15:35 AM Veronica Andreo wrote:
> Hi list
> 
> I have a vector, which are remotely sensed chlorophyll values for a 
certain
> pixel in 11 years, then i have a flag or label vector (t_max) that 
consists
> of 0 and 1, which tells me where i have the annual peak (one per 
year, then
> eleven 1 and all the rest are 0).
> 
> I'm interested in extracting the date in which the bloom beggins. For 
that
> matter I'm using a threshold of 5% above the median of the whole 
series.
> 
> I need to go over the chlorophyll data vector (cla) and every time I 
find
> the yearly maximum value (where t_max == 1), i need to go 
backwards in cla
> as many time steps as needed untill I cross the threshold (5% over 
the
> median of cl)... once i crossed that threshold, i need to be sure that 
at
> least for 2 time steps, the value of cla keeps going down (it's below 
the
> threshold)... once that condition is met, i want to get the first 
position
> where cla is higher than that threshold...
> that would be the phytoplankton bloom starting date...
> 
> This is more or less what i was trying, but it does not work and i'm
> stuck... can you please help me??
> 
> cla<- ts of 506 values (every 46 values, i have one year)
> t_max<- vector of 0 & 1; 1 where the yearly cla max occurs (11 
ones & all
> rest 0)
> threshold = 0.05*median(cla)+median(cla)
> 
> for ( t in 1:length(t_max) ) {
>   if ( t_max[t] == 0 )
>   {
> next;
>   } else {
>     if ( (cla[t-1] < threshold) && (cla[t-2)] < threshold) )
> {
>     return which(cla[t])
>     }
>   }
> }
> 
Hi Vero,
I think this does what you want, even though it doesn't do it the way 
you describe above.

cla<-runif(506)+
 rep(c(seq(1,3,length.out=23),seq(3,1,length.out=23)),11)
obs_time<-paste(rep(1:11,each=46),rep(1:46,11),sep="_")
t_max<-rep(0,506)
t_begin<-rep("",11)
for(year in 1:11) {
 threshold<-median(cla[(year-1) * 46 + 1:46]) * 1.05
 max_pos<-(year-1) * 46 + which.max(cla[(year-1)*46+1:46])
 t_max[max_pos]<-1
 threshpos<-(year-1) * 46 + 1
 while(cla[threshpos] < threshold) threshpos <- threshpos + 1
 cat(threshold,cla[threshpos],max_pos,cla[max_pos],"\n")
 t_begin[year]<-obs_time[threshpos]
}

Jim



More information about the R-help mailing list