[R-sig-Geo] gridded time series analysis

[Robert J. Hijmans]

Dear Advait,

> What I dont understand is:
> 1) How can each slope be the same (values of x1)

Apparently the values in these cells are all the same (not that you
only use the first 8 time steps, not all 120)

Here is an example with random values in which the slopes are all different

b <- brick(ncol=4, nrow=4, nl=120)
b <- setValues(b, matrix(runif(prod(dim(b))), ncol=120))
fun = function(x) { if(is.na(x[1])) { NA } else { lm(x[1:60] ~
x[61:120])$coefficients[2] } }
x <- calc(b, fun)
plot(x)
as.matrix(x)


> 2) The raster layers are of dimensions 4 x 4 x 120
> but dim(getValues(rcmraster_crop))
> [1]  16 120. Why is that? This must have implications for the way that the
> regressions for the brick are handled wouldnt it?
No

> Or is it specific to the
> way "getValues" handles things?

Yes. rows are cells (4*4=16), columns are layers.

> Also, does anyone have any pointers of how to accomplish the NA testing and
See an example above, assuming that if the first time step is NA, they
all are. You could also do things like
if (sum(is.na(x)) > 60)
depends on your context.


Best, Robert

On Mon, Dec 6, 2010 at 8:30 AM, Advait Godbole <advaitgodbole at gmail.com> wrote:
> Dear Members,
> To recap, I initially had trouble reading in two 3D netCDFs in to R using
> the raster "brick" function. Each brick is a 356 x 507 x 120 (lat, lon,
> time) raster. With that resolved (thanks Robert!) I have the following code
> to perform regressions between the  two time series (at each grid-cell)
> # read the netCDFs
> rcmraster <- brick("regcm.cumul.hdd.10k.96_05.fixed.latlonvars.nc", varname
> = "hdd")
> emitraster <-
> brick("flip.percap.vulcan.10k.dp.build.381.RES.EIAmy.96_05.this.latest.nc",varname="emitpercap")
> # crop subset of above rasters to test
> lim <- extent(rcmraster,180,185,300,305)  # createss a 4 x 4 x 120 brick
> rcmraster_crop <- crop(rcmraster,lim)
> emitraster_crop <- crop(emitraster,lim)
> #regressions
> s_crop <- stack(emitraster_crop,rcmraster_crop)
> func <- function(x) { lm(x[1:4] ~ x[5:8])$coefficients[2] }   #regression
> for emit* VS rcm*
> x1 <- calc(s_crop, func)
> While the above regressions run, I have to ultimately do the same  for the
>  entire domain (356 x 507 cells) and therefore need to check for NAs at all
> time steps for each grid cell for both the bricks, and then have lm run. To
> that  end,, I ran the following diagnostics on the cropped dataset:
>> dim(rcmraster_crop)
> [1]   4   4 120
>> dim(emitraster_crop)
> [1]   4   4 120
>> dim(x1)
> [1] 4 4
>> getValues(x1)
>  [1] 3.968207 3.968207 3.968207 3.968207 3.968207 3.968207 3.968207 3.968207
> 3.968207 3.968207 3.968207
> [12] 3.968207 3.968207 3.968207 3.968207 3.968207
>
> What I dont understand is:
> 1) How can each slope be the same (values of x1)
> 2) The raster layers are of dimensions 4 x 4 x 120
> but dim(getValues(rcmraster_crop))
> [1]  16 120. Why is that? This must have implications for the way that the
> regressions for the brick are handled wouldnt it? Or is it specific to the
> way "getValues" handles things?
> Also, does anyone have any pointers of how to accomplish the NA testing and
> regressions ?
> Thanks,
> Advait
>
> On Sun, Dec 5, 2010 at 6:29 PM, Robert J. Hijmans <r.hijmans at gmail.com>
> wrote:
>>
>> Advait,
>> It is just an example. In this case the function only checks the first
>> time-slice for all cells. That often works because oftentimes either
>> all values are NA, or none. Also, I think lm works when not all values
>> are NA. If your data is different, you may need a different solution.
>> The best thing might be to try it out with a small dataset (you can
>> use crop to create one).
>> Robert
>>
>> On Sun, Dec 5, 2010 at 3:14 PM, Advait Godbole <advaitgodbole at gmail.com>
>> wrote:
>> > Earlier in the threade you had mentioned that for two 120 layer bricks,
>> > this
>> > is what I should do:
>> > # regression of values in one brick (or stack) with another (Jacob's
>> > suggestion)
>> > s <- stack(s1, s2)
>> > # s1 and s2 have 120 layers
>> > fun <- function(x) { lm(x[1:120] ~ x[121:240])$coefficients[2] }
>> > x1 <- calc(s, fun)
>> > So if I have to check for NAs for the multi layered raster, I would need
>> > to
>> > do something like
>> > fun <- function (x) { if(is.na(x[1:120])) {NA} else
>> > {lm(x[1:120]~x[121:240])$coefficients[2]}}
>> > calc(s,fun=fun)
>> > That does not seem intuitively right to me although I am not very
>> > conversant
>> > with how is.na handles 3D variables to go much further. Also, should I
>> > not
>> > check the other brick for NAs too?
>> > Advait
>> > On Sun, Dec 5, 2010 at 5:29 PM, Robert J. Hijmans <r.hijmans at gmail.com>
>> > wrote:
>> >>
>> >> >If so, should I use a for loop instead?
>> >> No, that is the point of the calc function, that you do not do that.
>> >> Internally, calc uses apply as in
>> >> apply(values, 1, fun)
>> >>
>> >>
>> >> On Sat, Dec 4, 2010 at 1:50 PM, Advait Godbole
>> >> <advaitgodbole at gmail.com>
>> >> wrote:
>> >> > Roger, that fix worked, its there in the latest version. Thanks much.
>> >> > I am running into the problem that Martin mentioned earlier in the
>> >> > thread:
>> >> > Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok,
>> >> > ...) :
>> >> >  0 (non-NA) cases
>> >> >
>> >> > For which you suggested this:
>> >> > #you could catch this in your function and do something like:
>> >> >
>> >> >> fun=function(x) { if (is.na(x[1])){ NA } else {
>> >> >> lm(x[1]~x[2])$coefficients[2] }}
>> >> >> fun(c(NA,NA))
>> >> > [1] NA
>> >> >
>> >> > and then:
>> >> > calc(s, fun=fun)
>> >> >
>> >> > What I am confused is that my stack "s" that I will pass via calc is
>> >> > of
>> >> > dimensions 356 x 507 x 240, and therefore is.na will act on the
>> >> > entire
>> >> > array
>> >> > and therefore still fail at those points where it is not NA. Is that
>> >> > right?
>> >> > If so, should I use a for loop instead?
>> >> > Thanks,
>> >> > Advait
>> >> > On Fri, Dec 3, 2010 at 2:58 PM, Robert J. Hijmans
>> >> > <r.hijmans at gmail.com>
>> >> > wrote:
>> >> >>
>> >> >> It may have been fixed. Please try this again in raster 1.7-8 (now
>> >> >> on
>> >> >> CRAN for windows and linux, mac in a day or so). Robert
>> >> >>
>> >> >> On Fri, Dec 3, 2010 at 11:54 AM, steven mosher
>> >> >> <moshersteven at gmail.com>
>> >> >> wrote:
>> >> >> > Well,
>> >> >> >   I need some more information.
>> >> >> >  before you create the stack do this and report the results
>> >> >> > test1 <- layerNames(rcmraster)
>> >> >> > test2 <- layerNames(emitraster)
>> >> >> > str(test1)
>> >> >> > str(test2)
>> >> >> > Its a bit of a hack but  you can try to copy the layerNames out of
>> >> >> > the
>> >> >> > invidual stacks.
>> >> >> > Then set them to    ""
>> >> >> > then create one big stack
>> >> >> > then set the layerNames.
>> >> >> > ( I recall having to do that once in the past after a
>> >> >> > stackApply())
>> >> >> > I know in the past, I've had some similar issues, but I thought
>> >> >> > that
>> >> >> > they
>> >> >> >  had been fixed.
>> >> >> > I'll try to take a look at the code  in a bit, but if there is a
>> >> >> > bug,
>> >> >> > either
>> >> >> > robert or Jacob would have to fix it.
>> >> >> >
>> >> >> > On Fri, Dec 3, 2010 at 4:28 AM, Advait Godbole
>> >> >> > <advaitgodbole at gmail.com>
>> >> >> > wrote:
>> >> >> >>
>> >> >> >> With Robert being away until the 15th, can someone else help out
>> >> >> >> with
>> >> >> >> this
>> >> >> >> issue? To recap, I have read in two 507 x  356 x 120
>> >> >> >> (lon,lat,time)
>> >> >> >> netCDF
>> >> >> >> files into R using the "brick" function. I get the following
>> >> >> >> error
>> >> >> >> when
>> >> >> >> I
>> >> >> >> use "stack" to clump these two multi layer rasters together.
>> >> >> >> > s <- stack(rcmraster,emitraster)
>> >> >> >> Error in function (classes, fdef, mtable)  :
>> >> >> >>   unable to find an inherited method for function "trim", for
>> >> >> >> signature
>> >> >> >> "integer"
>> >> >> >> Robert suggested that it is because of the layerNames being
>> >> >> >> numeric,
>> >> >> >> annd
>> >> >> >> that layerNames() <- layerNames() would solve it but it did not
>> >> >> >> help,
>> >> >> >> even
>> >> >> >> though layerNames() is now of type character.
>> >> >> >> Thanks much!
>> >> >> >> adv
>> >> >> >> On Fri, Dec 3, 2010 at 7:13 AM, Advait Godbole
>> >> >> >> <advaitgodbole at gmail.com>
>> >> >> >> wrote:
>> >> >> >>>
>> >> >> >>> the setMinMax() function worked. Thanks Steve and Robert..
>> >> >> >>> the work around did not work - I still get the trim() error..
>> >> >> >>> adv
>> >> >> >>>
>> >> >> >>> On Thu, Dec 2, 2010 at 1:43 PM, Robert J. Hijmans
>> >> >> >>> <r.hijmans at gmail.com>
>> >> >> >>> wrote:
>> >> >> >>>>
>> >> >> >>>> This is the work around:
>> >> >> >>>>
>> >> >> >>>> layerNames(rcmraster) <- layerNames(rcmraster)
>> >> >> >>>>
>> >> >> >>>>
>> >> >> >>>> On Thu, Dec 2, 2010 at 10:43 AM, Robert J. Hijmans
>> >> >> >>>> <r.hijmans at gmail.com>
>> >> >> >>>> wrote:
>> >> >> >>>> > The question marks indicate that the min and max values are
>> >> >> >>>> > not
>> >> >> >>>> > known
>> >> >> >>>> > to the object (because the file type does not store them, and
>> >> >> >>>> > it
>> >> >> >>>> > is
>> >> >> >>>> > too 'expensive' to determine them on the fly). setMinMax if
>> >> >> >>>> > you
>> >> >> >>>> > really
>> >> >> >>>> > need to know.
>> >> >> >>>> > The error on trim seems to be a problem with layerNames being
>> >> >> >>>> > numeric.
>> >> >> >>>> > I'll fix that. This may be a work-around
>> >> >> >>>> >
>> >> >> >>>> >
>> >> >> >>>> >
>> >> >> >>>> >
>> >> >> >>>> >
>> >> >> >>>> >
>> >> >> >>>> > On Thu, Dec 2, 2010 at 10:23 AM, steven mosher
>> >> >> >>>> > <moshersteven at gmail.com> wrote:
>> >> >> >>>> >>  see the setMinMax()
>> >> >> >>>> >>
>> >> >> >>>> >>  for the error on trim() I'm not sure.
>> >> >> >>>> >>
>> >> >> >>>> >> On Thu, Dec 2, 2010 at 5:16 AM, Advait Godbole
>> >> >> >>>> >> <advaitgodbole at gmail.com>
>> >> >> >>>> >> wrote:
>> >> >> >>>> >>>
>> >> >> >>>> >>> That worked perfectly. Thanks for pointing this out. I
>> >> >> >>>> >>> changed
>> >> >> >>>> >>> the
>> >> >> >>>> >>> way my
>> >> >> >>>> >>> variables were being written to the netCDF and now "brick"
>> >> >> >>>> >>> works.
>> >> >> >>>> >>> However,
>> >> >> >>>> >>> I
>> >> >> >>>> >>> am getting an error while using the stack function, like
>> >> >> >>>> >>> so:
>> >> >> >>>> >>>
>> >> >> >>>> >>> "Error in function (classes, fdef, mtable)  :
>> >> >> >>>> >>>  unable to find an inherited method for function "trim",
>> >> >> >>>> >>> for
>> >> >> >>>> >>> signature
>> >> >> >>>> >>> "integer"
>> >> >> >>>> >>>
>> >> >> >>>> >>> The two rasters show the following:
>> >> >> >>>> >>> > rcmraster
>> >> >> >>>> >>> class       : RasterBrick
>> >> >> >>>> >>> filename    : regcm.cumul.hdd.10k.96_05.fixed.latlonvars.nc
>> >> >> >>>> >>> nlayers     : 120
>> >> >> >>>> >>> nrow        : 356
>> >> >> >>>> >>> ncol        : 507
>> >> >> >>>> >>> ncell       : 180492
>> >> >> >>>> >>> projection  : +proj=longlat +datum=WGS84 +ellps=WGS84
>> >> >> >>>> >>> +towgs84=0,0,0
>> >> >> >>>> >>> min value   : ? ? ? ? ? ? ? ? ? ?
>> >> >> >>>> >>> max value   : ? ? ? ? ? ? ? ? ? ?
>> >> >> >>>> >>> xmin        : -137.3312
>> >> >> >>>> >>> xmax        : -62.11259
>> >> >> >>>> >>> ymin        : 22.12693
>> >> >> >>>> >>> ymax        : 57.35954
>> >> >> >>>> >>> xres        : 0.1483602
>> >> >> >>>> >>> yres        : 0.09896801
>> >> >> >>>> >>>
>> >> >> >>>> >>> > emitraster
>> >> >> >>>> >>> class       : RasterBrick
>> >> >> >>>> >>> filename    :
>> >> >> >>>> >>>
>> >> >> >>>> >>>
>> >> >> >>>> >>>
>> >> >> >>>> >>> flip.percap.vulcan.10k.dp.build.381.RES.EIAmy.96_05.this.latest.nc
>> >> >> >>>> >>> nlayers     : 120
>> >> >> >>>> >>> nrow        : 356
>> >> >> >>>> >>> ncol        : 507
>> >> >> >>>> >>> ncell       : 180492
>> >> >> >>>> >>> projection  : +proj=longlat +datum=WGS84 +ellps=WGS84
>> >> >> >>>> >>> +towgs84=0,0,0
>> >> >> >>>> >>> min value   : ? ? ? ? ? ? ? ? ? ?
>> >> >> >>>> >>> max value   : ? ? ? ? ? ? ? ? ? ?
>> >> >> >>>> >>> xmin        : -137.3312
>> >> >> >>>> >>> xmax        : -62.11259
>> >> >> >>>> >>> ymin        : 22.12693
>> >> >> >>>> >>> ymax        : 57.35954
>> >> >> >>>> >>> xres        : 0.1483602
>> >> >> >>>> >>> yres        : 0.09896801
>> >> >> >>>> >>>
>> >> >> >>>> >>> A getValues call shows a lot of NAs. Now, my original file
>> >> >> >>>> >>> does
>> >> >> >>>> >>> have
>> >> >> >>>> >>> a
>> >> >> >>>> >>> large
>> >> >> >>>> >>> number of missing values so that is why this should happen,
>> >> >> >>>> >>> but I
>> >> >> >>>> >>> suspect
>> >> >> >>>> >>> something is going wrong since the min and max values are
>> >> >> >>>> >>> not
>> >> >> >>>> >>> being
>> >> >> >>>> >>> displayed for either of these rasters. Is that the source
>> >> >> >>>> >>> of
>> >> >> >>>> >>> the
>> >> >> >>>> >>> problem
>> >> >> >>>> >>> in
>> >> >> >>>> >>> the "stack" function?
>> >> >> >>>> >>>
>> >> >> >>>> >>> Thanks,
>> >> >> >>>> >>>
>> >> >> >>>> >>> Advait
>> >> >> >>>> >>>
>> >> >> >>>> >>> On Tue, Nov 30, 2010 at 5:18 PM, Robert J. Hijmans
>> >> >> >>>> >>> <r.hijmans at gmail.com>wrote:
>> >> >> >>>> >>>
>> >> >> >>>> >>> > Advait,
>> >> >> >>>> >>> >
>> >> >> >>>> >>> > The is an error in the file, I think, see below:
>> >> >> >>>> >>> >
>> >> >> >>>> >>> > > library(ncdf)
>> >> >> >>>> >>> > > filename =
>> >> >> >>>> >>> > > "c:/downloads/regcm.cumul.hdd.10k.96_05.new.nc"
>> >> >> >>>> >>> > > nc <- open.ncdf(filename)
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > zvar <- 'hdd'
>> >> >> >>>> >>> > > nc$var[[zvar]]$dim[[1]]$len
>> >> >> >>>> >>> > [1] 507
>> >> >> >>>> >>> > > nc$var[[zvar]]$dim[[2]]$len
>> >> >> >>>> >>> > [1] 356
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > xx <- nc$var[[zvar]]$dim[[1]]$vals
>> >> >> >>>> >>> > > # should be a vector of 507 longitudes but:
>> >> >> >>>> >>> > > dim(xx)
>> >> >> >>>> >>> > [1] 507 356 120
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > # works anyway
>> >> >> >>>> >>> > > xrange <- c(min(xx), max(xx))
>> >> >> >>>> >>> > > xrange
>> >> >> >>>> >>> > [1] -137.25705  -62.18677
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > yy <- nc$var[[zvar]]$dim[[2]]$vals
>> >> >> >>>> >>> > > # should be a vector of 356 latitudes but:
>> >> >> >>>> >>> > > dim(yy)
>> >> >> >>>> >>> > [1] 507 356 120
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > # this fails
>> >> >> >>>> >>> > > yrange <- c(min(yy), max(yy))
>> >> >> >>>> >>> > > yrange
>> >> >> >>>> >>> > [1] NaN NaN
>> >> >> >>>> >>> > > # and that is the only reason why raster failed.
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > #OK
>> >> >> >>>> >>> > > yy[,,1][1:10,350:355]
>> >> >> >>>> >>> >          [,1]     [,2]     [,3]     [,4]     [,5]
>> >> >> >>>> >>> > [,6]
>> >> >> >>>> >>> >  [1,] 51.49148 51.57122 51.65092 51.73055 51.81014
>> >> >> >>>> >>> > 51.88967
>> >> >> >>>> >>> >  [2,] 51.52888 51.60866 51.68839 51.76807 51.84769
>> >> >> >>>> >>> > 51.92727
>> >> >> >>>> >>> >  [3,] 51.56616 51.64598 51.72575 51.80547 51.88513
>> >> >> >>>> >>> > 51.96474
>> >> >> >>>> >>> >  [4,] 51.60332 51.68318 51.76299 51.84275 51.92245
>> >> >> >>>> >>> > 52.00210
>> >> >> >>>> >>> >  [5,] 51.64036 51.72026 51.80011 51.87991 51.95965
>> >> >> >>>> >>> > 52.03934
>> >> >> >>>> >>> >  [6,] 51.67728 51.75722 51.83711 51.91695 51.99673
>> >> >> >>>> >>> > 52.07646
>> >> >> >>>> >>> >  [7,] 51.71409 51.79407 51.87399 51.95387 52.03369
>> >> >> >>>> >>> > 52.11345
>> >> >> >>>> >>> >  [8,] 51.75077 51.83079 51.91076 51.99067 52.07053
>> >> >> >>>> >>> > 52.15033
>> >> >> >>>> >>> >  [9,] 51.78734 51.86739 51.94740 52.02735 52.10725
>> >> >> >>>> >>> > 52.18709
>> >> >> >>>> >>> > [10,] 51.82378 51.90388 51.98392 52.06391 52.14385
>> >> >> >>>> >>> > 52.22373
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > #bad values here (and note the NaN):
>> >> >> >>>> >>> > > yy[,,120][1:10,350:355]
>> >> >> >>>> >>> >                [,1]           [,2]           [,3]
>> >> >> >>>> >>> > [,4]
>> >> >> >>>> >>> >     [,5]           [,6]
>> >> >> >>>> >>> >  [1,] -1.313847e-229 -4.895880e-272 -7.553765e-266
>> >> >> >>>> >>> > 1.293814e+58
>> >> >> >>>> >>> > 9.551441e-293  -3.291580e-36
>> >> >> >>>> >>> >  [2,] -5.321424e+272 -3.076997e-223  -1.515780e-33
>> >> >> >>>> >>> > -9.895836e+171
>> >> >> >>>> >>> > 3.115686e-149  -4.970539e+92
>> >> >> >>>> >>> >  [3,]  5.259277e-203   5.980335e-56  2.893449e-212
>> >> >> >>>> >>> >  2.974200e+274
>> >> >> >>>> >>> > 3.186547e-218  2.503689e+296
>> >> >> >>>> >>> >  [4,] -3.977129e+225 -3.696719e-283   1.284612e+42
>> >> >> >>>> >>> > 1.887923e+28
>> >> >> >>>> >>> > -2.322520e+268  5.005820e-108
>> >> >> >>>> >>> >  [5,]  -6.569744e+88   1.824506e-74   7.663356e+39
>> >> >> >>>> >>> >  1.608759e-201
>> >> >> >>>> >>> > -1.433761e+110 -8.804183e-211
>> >> >> >>>> >>> >  [6,] -4.237822e-162            NaN  7.405904e-275
>> >> >> >>>> >>> >  5.859877e-258
>> >> >> >>>> >>> > 1.339820e+110  1.631993e+229
>> >> >> >>>> >>> >  [7,]   5.603325e-64   4.188764e+37 -2.188887e+243
>> >> >> >>>> >>> > -4.657907e+280
>> >> >> >>>> >>> > 1.795268e-185  1.166057e-175
>> >> >> >>>> >>> >  [8,] -2.112804e-229  4.536625e+137  2.910322e-278
>> >> >> >>>> >>> >  4.407645e-274
>> >> >> >>>> >>> > 8.886205e-239 -1.913860e+238
>> >> >> >>>> >>> >  [9,]  -9.126350e+44 -1.866422e+177 -2.578975e+256
>> >> >> >>>> >>> >  3.275379e+208
>> >> >> >>>> >>> > -2.657457e+75  -1.258259e-90
>> >> >> >>>> >>> > [10,]  6.698325e+216 -2.958355e+137 -1.058748e-178
>> >> >> >>>> >>> > -1.776776e+215
>> >> >> >>>> >>> > -6.368807e+78 -5.337653e+299
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > # this would have worked
>> >> >> >>>> >>> > > yrange <- c(min(yy[,,1]), max(yy[,,1]))
>> >> >> >>>> >>> > > yrange
>> >> >> >>>> >>> > [1] 22.17641 57.31006
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > close.ncdf(nc)
>> >> >> >>>> >>> >
>> >> >> >>>> >>> >
>> >> >> >>>> >>> > Best, Robert
>> >> >> >>>> >>> >
>> >> >> >>>> >>> >
>> >> >> >>>> >>> > On Tue, Nov 30, 2010 at 10:37 AM, Advait Godbole
>> >> >> >>>> >>> > <advaitgodbole at gmail.com> wrote:
>> >> >> >>>> >>> > > Dear Robert,
>> >> >> >>>> >>> > > here is the link to the file:
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > http://dl.dropbox.com/u/4030944/regcm.cumul.hdd.10k.96_05.new.nc. It
>> >> >> >>>> >>> > > is
>> >> >> >>>> >>> > a
>> >> >> >>>> >>> > > large file ~400 MB. Let me know if dropbox is not
>> >> >> >>>> >>> > > convenient
>> >> >> >>>> >>> > > and
>> >> >> >>>> >>> > > you
>> >> >> >>>> >>> > would
>> >> >> >>>> >>> > > rather do ftp -  but our machines are down for
>> >> >> >>>> >>> > > maintenance
>> >> >> >>>> >>> > > and I
>> >> >> >>>> >>> > > may
>> >> >> >>>> >>> > > not
>> >> >> >>>> >>> > be
>> >> >> >>>> >>> > > immediately able to get access.
>> >> >> >>>> >>> > > I was going over the NCL code that processes the files
>> >> >> >>>> >>> > > but
>> >> >> >>>> >>> > > did
>> >> >> >>>> >>> > > not
>> >> >> >>>> >>> > > find
>> >> >> >>>> >>> > > anything glaringly amiss, so I will have to look at it
>> >> >> >>>> >>> > > in
>> >> >> >>>> >>> > > greater
>> >> >> >>>> >>> > > detail
>> >> >> >>>> >>> > > since it is clear that the number of entries in the
>> >> >> >>>> >>> > > variables
>> >> >> >>>> >>> > > is
>> >> >> >>>> >>> > > not
>> >> >> >>>> >>> > right
>> >> >> >>>> >>> > > as you've pointed out. There is a chance that I have
>> >> >> >>>> >>> > > not
>> >> >> >>>> >>> > > written
>> >> >> >>>> >>> > > out
>> >> >> >>>> >>> > > the
>> >> >> >>>> >>> > > co-ordinates properly.
>> >> >> >>>> >>> > > Thanks,
>> >> >> >>>> >>> > > Advait
>> >> >> >>>> >>> > > On Mon, Nov 29, 2010 at 1:14 PM, Robert J. Hijmans
>> >> >> >>>> >>> > > <r.hijmans at gmail.com>
>> >> >> >>>> >>> > > wrote:
>> >> >> >>>> >>> > >>
>> >> >> >>>> >>> > >> 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
>> >> >> >>>> >>> > >> >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> > > --
>> >> >> >>>> >>> > > advait godbole
>> >> >> >>>> >>> > >
>> >> >> >>>> >>> >
>> >> >> >>>> >>>
>> >> >> >>>> >>>
>> >> >> >>>> >>>
>> >> >> >>>> >>> --
>> >> >> >>>> >>> advait godbole
>> >> >> >>>> >>>
>> >> >> >>>> >>>        [[alternative HTML version deleted]]
>> >> >> >>>> >>>
>> >> >> >>>> >>> _______________________________________________
>> >> >> >>>> >>> R-sig-Geo mailing list
>> >> >> >>>> >>> R-sig-Geo at r-project.org
>> >> >> >>>> >>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> >> >> >>>> >>
>> >> >> >>>> >>
>> >> >> >>>> >
>> >> >> >>>
>> >> >> >>>
>> >> >> >>>
>> >> >> >>> --
>> >> >> >>> advait godbole
>> >> >> >>
>> >> >> >>
>> >> >> >>
>> >> >> >> --
>> >> >> >> advait godbole
>> >> >> >
>> >> >> >
>> >> >
>> >> >
>> >> >
>> >> > --
>> >> > advait godbole
>> >> >
>> >
>> >
>> >
>> > --
>> > advait godbole
>> >
>
>
>
> --
> advait godbole
>