Re: [R-sig-eco] how to create box plot of precipitation data from raster stack
My short advice is to first treat all the files more abstractly, build a
data frame of the file name, and the full path to the file, and the date
associated with each (by string parse on the name).
Then you can check the set for sensibility, proper date, correct order,
sensible step between each. The full column of file paths can be given to
raster stack, and you can use setZ to apply the time layer reference in a
consistent way.
At another level up you might write functions that refer to you data frame,
and return layer/s or summaries as a function of date/season etc.
It helps to tidy up all the file stuff first, so it's a load of issues to
put behind you, and it's robust to future data coming in.
Cheers, Mike
On Sat, Mar 4, 2017, 14:26 karsten <kars...@terragis.net> wrote:
> Hi All,
>
> I am a beginner with R and spatial and working on plotting precipitation
> data (from CHIRPS) using raster stacks for the purpose of agricultural site
> characterization. I succeeded to generate plots for many research locations
> using a csv file with x-y coordinates and generating a bar graph as a
> 'climate diagram' for a selected year of a raster stack (containing several
> hundred *.tif files) with my r script.
> Now I would like to add a box plot to the output of each diagram showing
> the
> current years monthly rain data against the median, max mix of the entire
> 35
> years time series as a backdrop (what I would like to add should look like
> this http://www.terragis.net/docs/other/timeseries.png ). I am really
> struggling on how to prepare or extract data from the tifs in order to
> summarize them by month of the 35 years.
> Currently I have one folder with 420+ tif flies with names in this format
> es_af.1981.01.tif
> es_af.1981.02.tif
> es_af.1981.03.tif
> es_af.1981.04.tif
> es_af.1981.07.tif
> es_af.1981.08.tif
> es_af.1981.09.tif
> es_af.1981.10.tif
> es_af.1981.12.tif
> es_af.1982.01.tif
> and so on for 35 years
> How would I approach (write my r script) in order to extract a data frame /
> raster stack in order to output box plot diagrams as in the link above ?
> In my existing script I was using these steps to get a list of all file
> names and then creating individual stacks for each month , but what are my
> next steps ?
>
> # Use pattern arg to return a wildcard to get a list of all tifs in dir
> alltiffs = list.files(getwd(), pattern="*\\.tif$", full.names=TRUE)
>
> # get main data - e.g year 2016 to plot
> climategrids = alltiffs[grepl("*2016.*\\.tif$", alltiffs)]
>
> ## create list of all data files for all 35 years for each month
> ## then create stack for each month
> month_1_climategrids = alltiffs[grepl("*\\.01\\.tif$", alltiffs)]
> month_1 <- stack(month_1_climategrids, quick=TRUE)
> month_2_climategrids = alltiffs[grepl("*\\.02\\.tif$", alltiffs)]
> month_2 <- stack(month_2_climategrids, quick=TRUE)
>
> and so on
>
> Thanks
> Karsten Vennemann
> <http://www.terragis.net> www.terragis.net
>
>
> [[alternative HTML version deleted]]
>
> ___
> R-sig-ecology mailing list
> R-sig-ecology@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
>
--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia
[[alternative HTML version deleted]]
___
R-sig-ecology mailing list
R-sig-ecology@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
Re: [R-sig-eco] Extract values from 'tripGrid' output for Mantel test?
Hello, tripGrid returns a SpatialGridDataFrame (package sp) and so from that step you only need to know about that class. I'll use a data set in sp to create an object like yours: library(sp) data(meuse.grid) coordinates(meuse.grid) - ~x+y gridded(meuse.grid) - TRUE fullgrid(meuse.grid) - TRUE class(meuse.grid) tpc.for - meuse.grid class(tpc.for) [1] SpatialGridDataFrame attr(,package) [1] sp Notice that this object can contain multiple attributes, and this example data does: names(tpc.for) [1] part.a part.b dist soil ffreq but tripGrid just returns a grid with one attribute named z. You can index the layers by name or number e.g. tpc.for[1] or tpc.for[part.a], and plotting will assume just using the first. So to get data as an R matrix, you can use either of these, keeping in mind the different orientation conventions. I prefer the image() based one simply because I'm used to how it works, but they are just the transpose of each other: ff1 - as.matrix(tpc.for[1]) ff2 - as.image.SpatialGridDataFrame(tpc.for[1])$z (as.image... uses the same x,y,z name convention as image() which is partly why tripGrid unhelpfully uses that as well). Compare these which show that one is the column-reversed version of the other, and that image()'s orientation is not quite obvious: par(mfrow = c(2,2)) image(ff1, main = ff1 - as.matrix(x)) image(ff2, main = ff2 - as.image.SpatialGridDataFrame(x)$z) image(ff2[,ncol(ff2):1], main = ff2[,ncol(ff2):1]) I hope that helps, it probably does not matter which one you use but it's important to know the difference. Cheers, Mike. On Thu, Oct 4, 2012 at 9:51 PM, Kylie Scales klscale...@gmail.com wrote: Hi, I am testing levels of association between an environmental variable, extracted from a raster, and a time spent per unit area grid, created using the 'trip' package to determine whether the habtiat variable influences the amount of time spent in that grid cell by an animal. Both grids have been set to the same resolution. I would like to use a Mantel test to assess levels of 2D correlation. I have created a dissimilarity matrix using values extracted from the raster and now need a comparable matrix describing the 'time spent' data. These 'time spent' data are arranged as a slot in a SpatialGridDataFrame, known here as tpc.for@data$z. When I try to extract these data as a matrix: ff - as.matrix(tpc.for@data$z) the resultant object is only one-dimensional and then when I try to apply the function 'distance' from package 'ecodist', to calculate dissimilarity matrix, I get the following error: gg - distance(ff, bray-curtis) Error in numeric(N * N * P) : vector size cannot be NA In addition: Warning message: In N * N : NAs produced by integer overflow Any ideas on how I can extract the time spent values from the 'trip' object as a matrix, the dimensions of which will match the raster dimensions? Using 'image' plots the values at the same resolution. Any help would be very much appreciated. Kylie Scales [[alternative HTML version deleted]] ___ R-sig-ecology mailing list R-sig-ecology@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-ecology -- Michael Sumner Hobart, Australia e-mail: mdsum...@gmail.com ___ R-sig-ecology mailing list R-sig-ecology@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
Re: [R-sig-eco] How to calculate geographical coordinates of sampling points following (by bearing and distance) a georeferenced one?
I cannot remember the names, but there would be support for this in the maptools or geosphere package I am sure On Thursday, April 5, 2012, Ivailo wrote: Der fellow R-users, I have several sampling sites and in each site the starting plot (i.e. the first one that have been sampled) is georeferenced with a GPS but the following ones are just described by their distance and bearing in relation to the previous one. Is there a quick way to calculate (by using some of the spatial R-packages) the geographic coordinates of the sapling plots at each location that have not been explicitly georeferenced? Cheers, Ivailo -- UBUNTU: a person is a person through other persons. ___ R-sig-ecology mailing list R-sig-ecology@r-project.org javascript:; https://stat.ethz.ch/mailman/listinfo/r-sig-ecology -- Michael Sumner Institute for Marine and Antarctic Studies, University of Tasmania Hobart, Australia e-mail: mdsum...@gmail.com [[alternative HTML version deleted]] ___ R-sig-ecology mailing list R-sig-ecology@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
Re: [R-sig-eco] as.ltraj-adehabitatLT
The error message says it all - there are duplicated date/times in the da input for a given id. You should investigate these and either edit or remove them so that the trajectory makes sense. Something like this can show which ids/datetimes are duplicated, but without having your data it's just a guess (and untested): lapply(split(data.frame(da, data.sini$Name), data.sini$Name), function(x) x[is.duplicated(x),]) On Mon, Aug 8, 2011 at 6:58 AM, Natasa Djakovic natasa.djako...@umb.no wrote: Hi, I'm trying to create an object of class ltray to store animal movements, but keep getting an error message: Error in as.ltraj(xy = data.sini[, c(X, Y)], date = da, id = data.sini$Name) : non unique dates for a given burst. Does anyone have a suggestion how to solve this problem? Kind regards, Natasa ___ R-sig-ecology mailing list R-sig-ecology@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-ecology -- Michael Sumner Institute for Marine and Antarctic Studies, University of Tasmania Hobart, Australia e-mail: mdsum...@gmail.com ___ R-sig-ecology mailing list R-sig-ecology@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
Re: [R-sig-eco] How to calculate time spent on a world grid with circumpolar trips
Here's another way to do it, I've used a slightly different grid that is
only 360 cells across. This won't work if the trip segments also cross the 0
longitude, but I'll have a go at something in trip to do that too.
Cheers, Mike.
library(sp)
library(trip)
## read track data
a - read.table(file=79355.txt, h=T, sep=\t)
a$date - as.POSIXct(strptime(as.character(a$date),%d/%m/%Y %H:%M), GMT)
a$device_code - as.factor(a$device_code)
a - a[order(a$date), ]
## set up trip object
tr - a
tr$longitude[tr$longitude 0] - tr$longitude[tr$longitude 0] + 360
coordinates(tr) - ~longitude+latitude
tr - trip(tr, c(date, device_code))
## 0-360 grid
gt0360 - GridTopology(c(0, -70), c(1, 1), c(360, 57))
grd - tripGrid(tr, grid = gt0360)
## convert grid to -180/180 by indexing the x-coords and the grid
grd0 - as.image.SpatialGridDataFrame(grd)
subX - grd0$x 180
grd0$x[subX] - grd0$x[subX] - 360
grd0$x - sort(grd0$x)
grd0$z - rbind(grd0$z[subX,], grd0$z[!subX,])
## convert back to sp
grd - image2Grid(grd0)
On Fri, Oct 1, 2010 at 12:38 PM, Michael Sumner mdsum...@gmail.com wrote:
Thanks David, I'll investigate.
Also, I've been thinking there's probably no real reason not to allow
2-location trips. ;)
Cheers, Mike
2010/9/30 Pinaud David pin...@cebc.cnrs.fr:
Hi Mike!
I've fixed my problem by running this (unelegant) code, with the same
idea
as yours. I don't know if it's worth to include it in the trip package,
but
here is the code! Feel free to adapt it and use it if you want.
I needed to change a little bit the trip class definition to convert a
path
with two locations in a trip object. You set the limit at 3 locs when
defining a trip object, it makes sense that an interesting trip would
have
more than 3 locs, but some colleagues here are working with dispersal of
individuals, calculating distance and bearing with 2 locs (birth / death
or
birth/first reproduction place), so the package can be also usefull in
this
case.
Thanks
David
Le 09/09/2010 15:28, Michael Sumner a écrit :
Sorry to have missed this.
There's really no way to do this in longitude / latitude without a lot
of manual handling. You could split into two trips, one using the
eastern hemisphere, the other the western and use a shared [-180,180]
grid that you sum together.
What to do at the split coordinate depends on your actual data, there
are many options.
If you are still looking for an answer I could have a closer look at
your example code, but real data would be more helpful if you can
share it. There are a number of problems in the code you posted that
prevent it from being easily run.
Cheers, Mike.
On Fri, Aug 6, 2010 at 8:10 PM, Pinaud Davidpin...@cebc.cnrs.fr
wrote:
Dear all,
I want to calculate a grid of time spent per cell with albatross data.
The
goal is then to relate this to fisheries data set on a 1x1° grid.
I use the trip package with tripGrid(), but have some problems with
birds
crossing the date line.
Usually I use a projection to do that but here I need to keep the 1x1°
grid
in order to relate density tu some data (fisheries...).
I've tried with recenter() and nowrapSpatialLines(), without success.
Here a code:
library(trip)
library(maps) # to draw world maps
d- data.frame(x=c(seq(50, 180, 10), seq(-170, 50, 10)), y=rnorm(n=37,
m=-50, sd=3), tms=Sys.time()-(37:1*60*60*24), id=rep(1, 37)) # an
example
coordinates(d)- ~x+y
tr- trip(d, c(tms, id))
proj4string(at)- CRS(+proj=longlat +ellps=WGS84)
map('world')
plot(tr, add=T)
points(coordinates(tr), t=l)
# grid
grid.base- expand.grid(Long=seq(-180, 180, by=1), Lat=seq(-70, -15,
by=1))
# tab coords
coordinates(grid.base)- ~ Long + Lat
proj4string(grid.base)- CRS(+proj=longlat +ellps=WGS84)
gt- makeGridTopology(grid.base, cellsize = c(1, 1)) # def grid topo
# time spent per cell
tppc- tripGrid(tr, grid = gt)
tppc01- tppc # a copy for pres/abs
tpp...@data$z- ifelse(t...@data$z 0, 1, 0)
image(tppc, col=rev(heat.colors(500))) # seems ok
map('world', add=T)
points(coordinates(tr), t=l)
# but with presence/absence:
image(tppc01, col=c(white, green))
map('world', add=T)
points(coordinates(tr), t=l, col=blue)
Some ideas ?
Many thanks for your time.
David
--
***
David PINAUD
Ingénieur de Recherche Analyses spatiales
Centre d'Etudes Biologiques de Chizé - CNRS UPR1934
79360 Villiers-en-Bois, France poste 485
Tel: +33 (0)5.49.09.35.58
Fax: +33 (0)5.49.09.65.26
http://www.cebc.cnrs.fr/
***
__ Information from ESET Mail Security, version of virus
signature
database 5346 (20100806) __
The message was checked by ESET Mail Security.
http://www.eset.com
___
R-sig-ecology
