Jump to content


Photo

plot.lon.lat() in SGAT.

GeolocatorR SGAT FLightR plot.lon.lat()

Best Answer mdsumner , 13 October 2016 - 09:24 AM

Hi, it's been a while since I used SGAT directly so I need to refresh my memory a bit, but the output of Estelle is chains is samples from the posterior, and these are just arrays so it's not too hard to build up plots from them. I would try something like


## collapse out the parallel-structure

chain <- chainCollapse(fit$z)


xy <- apply(chain, 1:2, mean)

rn <- apply(chain, 1:2, sd)  ## sd may not be sensible, but other functions can be used here

plot(xy[,1], main = "longitude")

lines(xy[,1] - rn)

lines(xy[,2] + rn)


You'd want to get the range from rn and add it to the ylim first, something like

plot(xy[,1], main = "longitude", ylim = range(c(xy[,1] - rn, xy[,1] + rn))


and then do the same for latitude, on the same plot, etc. 


But, I tend to use Pimage to get more than than this rectangular "error box", and for that you can use any map projection (though this is all from memory). 


z <- Pimage(fit, type = "intermediate", proj = "laea")  

## laea Lambert Azimuthal Equal Area is a good starting guess, and by default this will build a local map projection



plot(z[1:20])  ## the time-weighted probability map for the first 20 time steps


z[5]  ## the  map for just the 5th time step


If you use type = "primary" you get the discrete estimates at the twilights ("x") rather than the smeared interval between them ("z"). 

With this you can employ any of the raster package tools to get proper contours of confidence intervals and so on, though we've never taken this to its full implementation. I'm happy to help if you contact me directly, or if you can share a fit object I can flesh out this here in more detail.  


Also I'd be happy to show how to get a plot like FlightR if you can show me one, it'd take a bit of work for me to do that myself. 


Cheers, Mike.  Go to the full post


  • Please log in to reply
3 replies to this topic

#1 AmelieRobertoCharron

AmelieRobertoCharron
  • Society Members
  • 9 posts
  • Winnipeg,
  • Canada

Posted 12 October 2016 - 02:37 PM

Hello! 

I am interested in plotting the results from the Estelle Model in SGAT as graphs of latitude and longitude over time while showing the error (similarly to the 'plot.lat.lon()' function in FLightR). Would anyone have any recommendations on how to do so? 

Any help would be appreciated! 

Thanks, 
Amélie



#2 ebridge

ebridge
  • Society Members
  • 27 posts

Posted 12 October 2016 - 05:23 PM

Hi Amélie

 

As a general approach you could find upper and lower quantiles for the posterior probability for longitude (stored in the output object fit$x) and latitude (fit$y) for each day and plot them along with the means or medians.

 

Are you looking for more of a "how to" guide?

 

Eli



#3 mdsumner

mdsumner
  • General Members
  • 3 posts
  • Tasmania,
  • Australia

Posted 13 October 2016 - 09:24 AM   Best Answer

Hi, it's been a while since I used SGAT directly so I need to refresh my memory a bit, but the output of Estelle is chains is samples from the posterior, and these are just arrays so it's not too hard to build up plots from them. I would try something like


## collapse out the parallel-structure

chain <- chainCollapse(fit$z)


xy <- apply(chain, 1:2, mean)

rn <- apply(chain, 1:2, sd)  ## sd may not be sensible, but other functions can be used here

plot(xy[,1], main = "longitude")

lines(xy[,1] - rn)

lines(xy[,2] + rn)


You'd want to get the range from rn and add it to the ylim first, something like

plot(xy[,1], main = "longitude", ylim = range(c(xy[,1] - rn, xy[,1] + rn))


and then do the same for latitude, on the same plot, etc. 


But, I tend to use Pimage to get more than than this rectangular "error box", and for that you can use any map projection (though this is all from memory). 


z <- Pimage(fit, type = "intermediate", proj = "laea")  

## laea Lambert Azimuthal Equal Area is a good starting guess, and by default this will build a local map projection



plot(z[1:20])  ## the time-weighted probability map for the first 20 time steps


z[5]  ## the  map for just the 5th time step


If you use type = "primary" you get the discrete estimates at the twilights ("x") rather than the smeared interval between them ("z"). 

With this you can employ any of the raster package tools to get proper contours of confidence intervals and so on, though we've never taken this to its full implementation. I'm happy to help if you contact me directly, or if you can share a fit object I can flesh out this here in more detail.  


Also I'd be happy to show how to get a plot like FlightR if you can show me one, it'd take a bit of work for me to do that myself. 


Cheers, Mike. 

#4 AmelieRobertoCharron

AmelieRobertoCharron
  • Society Members
  • 9 posts
  • Winnipeg,
  • Canada

Posted 19 October 2016 - 09:33 AM

Thank you very much for the help! The code you provided worked perfectly!







Also tagged with one or more of these keywords: GeolocatorR, SGAT, FLightR, plot.lon.lat()

0 user(s) are reading this topic

0 members(s), 0 guests(s) and 0 anonymous member(s)