ctmm-initiative / ctmm

Continuous-Time Movement Modeling. Functions for identifying, fitting, and applying continuous-space, continuous-time stochastic movement models to animal tracking data.
http://biology.umd.edu/movement.html
43 stars 10 forks source link

`ctmm::simulate`, how to backtransform coordinates? #28

Closed jmsigner closed 5 years ago

jmsigner commented 5 years ago

I am having some troubles understanding how to back transform coordinates to the original coordinates system after simulating a path. I tried to illustrate the problem with a small repex (below).

Any hints or pointers to what I am doing wrong would be greatly be appreciated.

library(lubridate)
library(move)
library(ctmm)

set.seed(123)

# Generate some dummy data
dat <- data.frame(x = cumsum(rnorm(100)),
                  y = cumsum(rnorm(100)),
                  t = ymd("2018-01-01") + hours(0:99))

# Create a telemetry object via move
m <- move(x = dat$x, y = dat$y,
                time = dat$t, proj = sp::CRS("+init=epsg:4326"),
                data = data.frame(individual.local.identifier = rep("1", nrow(dat))))

c1 <- as.telemetry(m)

# Fit a ou model
g <- ctmm.guess(c1, interactive = FALSE)
ou <- ctmm.fit(c1, ctmm(tau = g$tau[1]))

# Simulate from this model
t <- as.numeric(dat$t)
s1 <- simulate(ou, data = c1, t = t)

# Coordinates are very differnt
plot(s1)

# from
plot(dat[, 1:2])

# Coerce to spatial points doesent help
p1 <- SpatialPointsDataFrame(s1)
sp::plot(p1, axes = TRUE)

# from
plot(dat[, 1:2])
jmsigner commented 5 years ago

I am using:

ctmm              * 0.5.3      2018-10-09 Github (ctmm-initiative/ctmm@46ee5ea)
chfleming commented 5 years ago

I can see the issue. The move package is not strict about what projections it will allow. When as.telemetry() sees "+init=epsg:4326 +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0" from projection(m), it infers this projection not to be in meters but degrees and generates an automated projection (locally flat, in meters) so that analysis can proceed in a sensible way. You can see that it is in a different projection from projection(c1). p1 then inherits its projection from c1, after running through SpatialPointsDataFrame.telemetry(). I will add a warning when this happens in as.telemetry().

I'm not sure how to reproject c1 from within sp. I tried projection()<- and for me that just overwrote the projection string without actually reprojecting the data. I am not very familiar with sp, but provide the export facilities on request.

jmsigner commented 5 years ago

If you use sp::spTransform the the CRS is actually transformed not just overwritten (as it is the case with projection).

It works well, if I omit the move package and use 'relatively small' coord values.

library(ctmm)
set.seed(123)

# Generate some dummy data: works well
dat <- data.frame(lon = cumsum(rnorm(100)),
                  lat = cumsum(rnorm(100)),
                  timestamp = ymd("2018-01-01") + hours(0:99))

xx <- as.telemetry(dat)
slot(xx, ".Data")

# Generate some dummy data: coords are Inf, which causes troubles later
# I assume as.telemetry expects geographical coordinates
dat <- data.frame(lon = 1e3 + cumsum(rnorm(100)),
                  lat = 1e3 + cumsum(rnorm(100)),
                  timestamp = ymd("2018-01-01") + hours(0:99))

xx <- as.telemetry(dat)
slot(xx, ".Data")

# Generate some dummy data: coords are Inf, which causes troubles later 
# But even if I provide a different CRS, it appears not be recognized.
s <- SpatialPoints(cbind(10, 52), CRS("+init=epsg:4326"))
s <- coordinates(spTransform(s, CRS("+init=epsg:3035")))
dat <- data.frame(lon = s[1, 1] + cumsum(rnorm(100)),
                  lat = s[1, 2] + cumsum(rnorm(100)),
                  timestamp = ymd("2018-01-01") + hours(0:99))

xx <- as.telemetry(dat, projection = CRS("+init=epsg:3035"))
slot(xx, ".Data")
chfleming commented 5 years ago

as.telemetry() assumes Movebank format, so, yes, geographical coordinates.

chfleming commented 5 years ago

The latest update should have fixed the projection to be preserved if valid (it wasn't before) and a warning issued otherwise.

briana-abrahms commented 5 years ago

Hi Chris,

I'm still confused on this issue. I have my initial coordinates in long/lat. When I simulate on this the simulated telemetry object seems to report back the projection system I specified, but the data in the simulated object is still in flat meters. Example

s <- SpatialPoints(cbind(-80, 30), CRS("+init=epsg:4326"))
dat <- data.frame(lon = coordinates(s)[1,1] + cumsum(rnorm(100)),
                  lat = coordinates(s)[1,2] + cumsum(rnorm(100)),
                  timestamp = ymd("2018-01-01") + hours(0:99))

t <- as.telemetry(dat, projection = "+init=epsg:4326")
slot(t, ".Data") #6 slots: YYYY-MM-DD HH:MM:SS, lat, lon, t, x, y
plot(t) #flat meters
projection(t) #"+init=epsg:4326"

# Fit ou model
g <- ctmm.guess(t, interactive = FALSE)
ou <- ctmm.fit(t, ctmm(tau = g$tau[1]))

# Simulate from this model
s1 <- simulate(ou, data = t, t = dat$dt)
slot(s1, ".Data") #3 slots: t, x, y
plot(s1) #flat meters
projection(s1) "+init=epsg:4326"

How would one back-transform t,x,y to yyyy-mm-dd hh:mm:ss, lon, lat?

Thanks, Briana

chfleming commented 5 years ago

Hi Briana,

The first thing I see is that I don't have a sufficient check on the projection argument. projection = "+init=epsg:4326" is long-lat and so as.telemetry should throw an error there. (ctmm calculations need a locally flat projection in meters.) I can fix that by running the projection string through CRS to put it into canonical form before checking for longlat/latlong.

The second thing is that I suppose you want a convenience function to populate simulate/predict telemetry objects with proper timestamps and long-lat? I can do that. I just need a day or two to get to it.

Thanks, Chris

briana-abrahms commented 5 years ago

Hi Chris,

Thanks for your quick response. Yes, that would be great re: populating telemetry objects with the full timestamps and long-lat.

Best, Briana

On Sat, Nov 10, 2018 at 10:10 AM Christen H. Fleming < notifications@github.com> wrote:

Hi Briana,

The first thing I see is that I don't have a sufficient check on the projection argument. projection = "+init=epsg:4326" is long-lat and so as.telemetry should throw an error there. (ctmm calculations need a locally flat projection in meters.) I can fix that by running the projection string through CRS to put it into canonical form before checking for longlat/latlong.

The second thing is that I suppose you want a convenience function to populate simulate/predict telemetry objects with proper timestamps and long-lat? I can do that. I just need a day or two to get to it.

Thanks, Chris

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ctmm-initiative/ctmm/issues/28#issuecomment-437604173, or mute the thread https://github.com/notifications/unsubscribe-auth/AWKptl6CR-FCPI_g9X1p3QZYosIoU_l9ks5utxZ4gaJpZM4XafV6 .

chfleming commented 5 years ago

Hi Briana,

In the development branch, simulate and predict now have a complete argument that can be used to calculate the timestamps and geographic coordinates.

https://ctmm-initiative.github.io/ctmm/reference/simulate.ctmm.html

briana-abrahms commented 5 years ago

Awesome! Thank you so much!

On Tue, Nov 13, 2018 at 2:40 PM Christen H. Fleming < notifications@github.com> wrote:

Hi Briana,

In the development branch, simulate and predict now have a complete argument that can be used to calculate the timestamps and geographic coordinates.

https://ctmm-initiative.github.io/ctmm/reference/simulate.ctmm.html

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ctmm-initiative/ctmm/issues/28#issuecomment-438466127, or mute the thread https://github.com/notifications/unsubscribe-auth/AWKptuzeAGT0M4GqUlrNWxOMCpRqC97Oks5uu0phgaJpZM4XafV6 .