[R-sig-Geo] gridded time series analysis

Robert J. Hijmans r.hijmans at gmail.com
Mon Nov 29 19:14:53 CET 2010


Advait,
I do not think I can help you further unless you send me your file;
but this suggests that it is either not following the CF standards, or
something is wrong with the file, or that there is something I have
not understood about these files.
You said you have 507 x 356 x 120 values. I would expect
length(frcm$dim[[1]]$vals) == 507 and length(frcm$dim[[2]]$vals) ==
356
Robert

On Mon, Nov 29, 2010 at 6:11 AM, Advait Godbole <advaitgodbole at gmail.com> wrote:
> sorry to disrupt the flow of the conversation, but after doing what Robert
> suggested, this is what was returned:
>> length(frcm$dim[[1]]$vals)
> [1] 21659040
>> length(frcm$dim[[2]]$vals)
> [1] 21659040
> the other two functions return the latitude and longitude arrays not
> populated with NAs/missing values. I havent reproduced them here because of
> their number.
> Regards,
> advait
> On Sat, Nov 27, 2010 at 3:48 PM, Robert J. Hijmans <r.hijmans at gmail.com>
> wrote:
>>
>> Advait,
>>
>> That suggests that you have missing values in the lon or lat
>> dimensions. Is that true? That would be awkward, but perhaps something
>> else is going on. Can you send me your file? If not can you send the
>> results of this:
>>
>> filename <- "regcm.cumul.hdd.10k.96_05.new.nc"
>> nc <- open.ncdf(filename)
>> length(nc$dim[[1]]$vals)
>> length(nc$dim[[2]]$vals)
>> nc$dim[[1]]$vals
>> nc$dim[[2]]$vals
>> close.ncdf(nc)
>>
>> Thanks, Robert
>>
>> On Sat, Nov 27, 2010 at 8:36 AM, Advait Godbole <advaitgodbole at gmail.com>
>> wrote:
>> > Thanks everyone for these leads..However, I am running into problems
>> > with
>> > data handling. I used ncdf to extract the two variables (507 x 356 x 120
>> > -
>> > lon,lat,time) and converted them into array objects. The suggestions
>> > provided all utilize the raster package, but using "brick" on the file
>> > returns this error: "rcmraster <-
>> > brick("regcm.cumul.hdd.10k.96_05.new.nc")
>> > Error in if (e at xmin > -400 & e at xmax < 400 & e at ymin > -90.1 & e at ymax <  :
>> > missing value where TRUE/FALSE needed"
>> >
>> > I have a fair number of NAs (all missing values in the original
>> > variables).
>> > Is there a way to set a flag for that? Or any other workaround?
>> > Thanks much,
>> > advait
>> > On Fri, Nov 26, 2010 at 6:25 PM, Robert J. Hijmans <r.hijmans at gmail.com>
>> > wrote:
>> >>
>> >> There are some difference in the behavior of 'calc' between functions,
>> >> because of attempts to accommodate different functions & intentions.
>> >> But in 'raster' 1.7-4 (available from R-Forge in ~ 24 hrs; and from
>> >> CRAN soon), the below will work:
>> >>
>> >> library(raster)
>> >> #creating some data
>> >> r <- raster(nrow=10, ncol=10)
>> >> s1 <- s2<- list()
>> >> for (i in 1:12) {
>> >>        s1[i] <- setValues(r, rnorm(ncell(r), i, 3) )
>> >>        s2[i] <- setValues(r, rnorm(ncell(r), i, 3) )
>> >> }
>> >> s1 <- stack(s1)
>> >> s2 <- stack(s2)
>> >>
>> >> # regression of values in one brick (or stack) with another (Jacob's
>> >> suggestion)
>> >> s <- stack(s1, s2)
>> >> # s1 and s2 have 12 layers
>> >> fun <- function(x) { lm(x[1:12] ~ x[13:24])$coefficients[2] }
>> >> x1 <- calc(s, fun)
>> >>
>> >> # regression of values in one brick (or stack) with 'time'
>> >> time <- 1:nlayers(s)
>> >> fun <- function(x) { lm(x ~ time)$coefficients[2] }
>> >> x2 <- calc(s, fun)
>> >>
>> >> # get multiple layers, e.g. the slope _and_ intercept
>> >> fun <- function(x) { lm(x ~ time)$coefficients }
>> >> x3 <- calc(s, fun)
>> >>
>> >> If this does not work in your version, you can use apply( ) as in what
>> >> I send earlier.
>> >>
>> >> Robert
>> >>
>> >> On Fri, Nov 26, 2010 at 2:56 PM, Robert J. Hijmans
>> >> <r.hijmans at gmail.com>
>> >> wrote:
>> >> > It seems that 'calc' does not like this (any more; another thing to
>> >> > fix) . If your rasters are not very large you can use apply, which
>> >> > makes it much faster:
>> >> >
>> >> > library(raster)
>> >> > #creating some data
>> >> > r <- raster(nrow=10, ncol=10)
>> >> > s <- list()
>> >> > for (i in 1:25) { s[i] <- setValues(r, rnorm(ncell(r), i, 3) ) }
>> >> > s <- stack(s)
>> >> >
>> >> > # regression
>> >> > time <- 1:nlayers(s)
>> >> > fun <- function(x) { lm(x ~ time)$coefficients[2] }
>> >> > r[] <- apply(getValues(s), 1, fun)
>> >> >
>> >> > Robert
>> >> >
>> >> >
>> >> >
>> >> > On Fri, Nov 26, 2010 at 2:51 PM, Jacob van Etten
>> >> > <jacobvanetten at yahoo.com> wrote:
>> >> >> you could try this approach (use calc whenever you can):
>> >> >>
>> >> >> (supposing your bricks have 12 layers)
>> >> >>
>> >> >> br3 <- stack(brick1, brick2)
>> >> >> lmS <- function(x) lm(x[1:12] ~ x[13:24)$coefficients[2]
>> >> >> r <- calc(br3, lmS)
>> >> >>
>> >> >> Jacob.
>> >> >>
>> >> >> --- On Fri, 26/11/10, steven mosher <moshersteven at gmail.com> wrote:
>> >> >>
>> >> >> From: steven mosher <moshersteven at gmail.com>
>> >> >> Subject: Re: [R-sig-Geo] gridded time series analysis
>> >> >> To: "Martin" <martin_brandt at gmx.net>
>> >> >> Cc: r-sig-geo at stat.math.ethz.ch
>> >> >> Date: Friday, 26 November, 2010, 23:33
>> >> >>
>> >> >> that's cool, I'm also interested in a similar problem but just with
>> >> >> one
>> >> >> brick ending up with a slope raster as the output. It may be
>> >> >> possible
>> >> >> with
>> >> >> stackApply(). have a look. or maybe robert will chime in
>> >> >>
>> >> >>
>> >> >>
>> >> >> On Fri, Nov 26, 2010 at 1:35 PM, Martin <martin_brandt at gmx.net>
>> >> >> wrote:
>> >> >>
>> >> >>>
>> >> >>> this is what I did to perform a regression between two bricks (each
>> >> >>> brick
>> >> >>> represent a time series):
>> >> >>>
>> >> >>> r <- raster(brick1)
>> >> >>> for (i in 1:ncell(r)) {
>> >> >>> r[i] = lm(as.ts(cellValues(brick1, i)) ~ as.ts(cellValues(brick2,
>> >> >>> i)))$coefficients[2]
>> >> >>> }
>> >> >>>
>> >> >>> The result will be a slope raster, but it really takes a lot of
>> >> >>> time,
>> >> >>> so
>> >> >>> maybe there is a better solution..
>> >> >>>
>> >> >>>
>> >> >>> --
>> >> >>> View this message in context:
>> >> >>>
>> >> >>>
>> >> >>> http://r-sig-geo.2731867.n2.nabble.com/gridded-time-series-analysis-tp5775651p5778472.html
>> >> >>> Sent from the R-sig-geo mailing list archive at Nabble.com.
>> >> >>>
>> >> >>> _______________________________________________
>> >> >>> R-sig-Geo mailing list
>> >> >>> R-sig-Geo at stat.math.ethz.ch
>> >> >>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >> >>>
>> >> >>
>> >> >>     [[alternative HTML version deleted]]
>> >> >>
>> >> >> _______________________________________________
>> >> >> R-sig-Geo mailing list
>> >> >> R-sig-Geo at stat.math.ethz.ch
>> >> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >> >>
>> >> >>
>> >> >>
>> >> >>
>> >> >>        [[alternative HTML version deleted]]
>> >> >>
>> >> >>
>> >> >> _______________________________________________
>> >> >> R-sig-Geo mailing list
>> >> >> R-sig-Geo at stat.math.ethz.ch
>> >> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >> >>
>> >> >>
>> >> >
>> >>
>> >> _______________________________________________
>> >> R-sig-Geo mailing list
>> >> R-sig-Geo at stat.math.ethz.ch
>> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >
>> >
>> >
>> > --
>> > advait godbole
>> >
>
>
>
> --
> advait godbole
>



More information about the R-sig-Geo mailing list