Re: [R-sig-Geo] Finding the Circumscribing Circle of a polygon
Thanks Michael for the answer. How come that I didn't find this solution on google?...When you use rseek and type circumscribing circle, you don't find any mention to the tripack package..Maybe it's because the fnction has been developed lately? the latest version of tripack being from february 2013 Anyway, it's nice to have a solution, now. What's the best source when you look for a function to accomplish something? 2013/5/16 Michael Sumner > Just FYI, there is tripack::circumcircle, here's an example with data > in the maptools package. > > library(maptools) > > data(wrld_simpl) > > poly <- wrld_simpl[wrld_simpl$NAME == "Australia", ] > > ## use coercion to get raw vertices (long winded, but avoids @ usage) > coords <- coordinates(as(as(poly, "SpatialLines"), "SpatialPoints")) > > > library(tripack) > cc<- circumcircle(coords[,1], coords[,2], plot = TRUE) > plot(poly, add = TRUE, col = "grey") > > > On Fri, May 17, 2013 at 7:39 AM, JLong wrote: > > Hey Mathieu and Greg, > > > > Great ideas, and thanks for the code Mathieu. > > I can't imagine a polygon shape where the 2nd step would ever by > necessary? > > Can you elaborate on why that might be. > > > > FYI, Mathieu I'm attempting to compute the same shape/linearity metric as > > you. > > > > Cheers, > > Jed > > > > > > > > > > - > > Jed Long > > PhD Candidate > > SPAR Lab, Dept. of Geography > > University of Victoria, Canada > > -- > > View this message in context: > http://r-sig-geo.2731867.n2.nabble.com/Finding-the-Circumscribing-Circle-of-a-polygon-tp7582372p7583598.html > > Sent from the R-sig-geo mailing list archive at Nabble.com. > > > > ___ > > R-sig-Geo mailing list > > R-sig-Geo@r-project.org > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > > > > -- > Michael Sumner > Hobart, Australia > e-mail: mdsum...@gmail.com > > ___ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Finding the Circumscribing Circle of a polygon
Just FYI, there is tripack::circumcircle, here's an example with data in the maptools package. library(maptools) data(wrld_simpl) poly <- wrld_simpl[wrld_simpl$NAME == "Australia", ] ## use coercion to get raw vertices (long winded, but avoids @ usage) coords <- coordinates(as(as(poly, "SpatialLines"), "SpatialPoints")) library(tripack) cc<- circumcircle(coords[,1], coords[,2], plot = TRUE) plot(poly, add = TRUE, col = "grey") On Fri, May 17, 2013 at 7:39 AM, JLong wrote: > Hey Mathieu and Greg, > > Great ideas, and thanks for the code Mathieu. > I can't imagine a polygon shape where the 2nd step would ever by necessary? > Can you elaborate on why that might be. > > FYI, Mathieu I'm attempting to compute the same shape/linearity metric as > you. > > Cheers, > Jed > > > > > - > Jed Long > PhD Candidate > SPAR Lab, Dept. of Geography > University of Victoria, Canada > -- > View this message in context: > http://r-sig-geo.2731867.n2.nabble.com/Finding-the-Circumscribing-Circle-of-a-polygon-tp7582372p7583598.html > Sent from the R-sig-geo mailing list archive at Nabble.com. > > ___ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo -- Michael Sumner Hobart, Australia e-mail: mdsum...@gmail.com ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Finding the Circumscribing Circle of a polygon
Hey Mathieu and Greg, Great ideas, and thanks for the code Mathieu. I can't imagine a polygon shape where the 2nd step would ever by necessary? Can you elaborate on why that might be. FYI, Mathieu I'm attempting to compute the same shape/linearity metric as you. Cheers, Jed - Jed Long PhD Candidate SPAR Lab, Dept. of Geography University of Victoria, Canada -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Finding-the-Circumscribing-Circle-of-a-polygon-tp7582372p7583598.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Spatially autocorrelated random effects in logistic regression.
Hi Seth,
Thanks for your input on this. Unfortunately it doesn't seem like it will
accept 2 grouping factors. In fact, it doesn't let you specify different groups
if your response is already grouped, i.e. as a matrix of positives and
negatives. For example
Y<-cbind(1:5,5:1)Long<-c(0,1,1,0,1)Lat<-c(0,0,1,2,2)Covar<-seq(10,50,10)Group=factor(rep("a",nrow(Y)))g<-factor(1:5)
All of these produce errors:
# run a model just with random=argument for householdmodel1 = glmmPQL(Y ~
Covar, random = ~ 1 | g, correlation=corExp(form=~Long+Lat),
family=binomial(logit))
model2 = glmmPQL(Y ~ Covar, random = ~ 1 | g,
correlation=corExp(form=~Long+Lat|g), family=binomial(logit))
model3 = glmmPQL(Y ~ Covar, random = ~ 1 | g,
correlation=corExp(form=~Long+Lat|Group), family=binomial(logit)) ## This is
the one I really want..!
# The only model that runs ismodel4 = glmmPQL(Y ~ Covar, random = ~ 1 |
Group, correlation=corExp(form=~Long+Lat), family=binomial(logit))
# This seems to be the same as model4, i.e. model4 is borrowing the group level
specified as random in the correlation model5 = glmmPQL(Y ~ Covar, random = ~
1 | Group, correlation=corExp(form=~Long+Lat | Group), family=binomial(logit))
# It doesn't seem to like the corMatrix command on glmmPQL
modelscorMatrix(model5) Error in UseMethod("corMatrix") : no applicable
method for 'corMatrix' applied to an object of class "c('glmmPQL', 'lme')"
I've done a fair amount of reading on the use of glmmPQL but nowhere can I find
anything on this specific problem..
Date: Thu, 16 May 2013 15:27:56 -0400
Subject: Re: [R-sig-Geo] Spatially autocorrelated random effects in logistic
regression.
From: sjmyers3...@gmail.com
To: hughsturr...@hotmail.com
Hi,
I forgot about the bit of "tricking" glmmPQL by adding a dummy variable that is
the same for each observation. Do you have SAS by any chance? This may be a
little easier in Proc Glimmix.
I believe you can have two grouping factors. You can have one grouping factor
where each household has its own grouping, so household one is "1" (as a
factor) and household 2 is "2" and there would be multiple observations within
each nesting at that grouping level. This is what you would specify in your
random = argument if you want to use a random intercept to account for
correlations within a single household. You can then specify another grouping
factor where all observations are a "1" or some other factor. Then you can use
the | symbol within the correlation= argument after the east + north to model
between household spatial autocorrelation. It seems that glmmPQL might be
"borrowing" the grouping factor for correlation= from the random=, which is odd
if true (I just skimmed an attached chapter, and this seems to be the case, but
perhaps you can force two different grouping factors one at the household and
one where the entire study area is a group).
Try specifying random as you had it with house as the group and then
correlation as
correlation=corExp(form=~Long+Lat | Group)
I believe the partial sill and nugget effect won't do anything if the
coordinates are exactly the same. I am not 100% certain though. It is down to
how it was coded. You can post this question to the appropriate mailing list
and maybe get the package authors to respond. I would try it though to see in
the following way.
run a model just with random=argument for householdrun a model just with
correlation=argument with all data in one big group (no nugget)run a model just
with correlation=argument with all data in one big group (WITH nugget)
run a model just with random=argument for household AND with
correlation=argument with all data in one big group (no nugget)
run a model just with random=argument for household AND with
correlation=argument with all data in one big group (WITH nugget)
If you save each model and look at the correlation matrices for each and
compare them, you can likely see what is going on.
Here's a book chapter that I read when I first ran glmmPQL. The relevant part
starts around page 301.
Seth
On Thu, May 16, 2013 at 2:06 PM, Hugh Sturrock wrote:
Yes that's what I thought. I need to specify a different group for each
household, however, as your pasted section states, if a grouping factor is
specified, the correlation structure is assumed to apply to observations WITHIN
the grouping level, not BETWEEN. I would like to do both, i.e. account for
correlation within groups and spatial autocorrelation between groups.
If I run the model with the response as a vector of 1's and 0's I have the
problem of cooincident coordinates. If I run the model with the response as a
matrix of positives and negatives, when I specify each row of the response as a
separate group, it gives me the errors I posted before. The only way seems to
be to force each row of the reponse (which is several observations of positives
and negatives at each site) to be a single 'observation', create one
Re: [R-sig-Geo] nearest neighbour using distance in metres
Ross,
The geodist function in the package GMT will compute distance in all sorts of
units (m, km etc.) from lat longs.
Hugh
> Date: Thu, 16 May 2013 14:16:03 +0100
> From: rossah...@googlemail.com
> To: r-sig-geo@r-project.org
> Subject: [R-sig-Geo] nearest neighbour using distance in metres
>
> Hi all,
>
> Following code computes Euclidian distance from each point to it nearest
> neighbour:
>
> library(spatstat)
> y <- rnorm(30, 55.65, .008)
> x <- rnorm(30, -1.82, .01)
> xy <- cbind(x, y)
> nndist(xy)
>
> Is it possible to use distance in metres, not Euclidian distances? Here are
> same points projected:
>
> xySpatial <- SpatialPoints(xy, proj4string= CRS("+proj=longlat +datum=WGS84
> +ellps=WGS84 +towgs84=0,0,0"))
> xySpatialTransformed <- spTransform(xySpatial, CRS("+proj=laea +lon_0=-1.92
> +lat_0=55.7 +datum=WGS84"))
>
> Thanks
> Ross
>
>
>
> [[alternative HTML version deleted]]
>
> ___
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[[alternative HTML version deleted]]
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] hand rolled spatial correlation in gnls and how to extract conditioned residuals?
I'm not getting any help from r-sig-mixed mailing list, so I thought I'd ask this here as it is related to spatial modeling. Questions: 1. Is there anything obviously wrong in concept or code with the following approach? 2. How do I extract residuals from models using correlation= argument so that the residuals are "adjusted" by the correlation among residuals that is specified? Background: I have a nonlinear model with a continuous response variable (0<=resp.>=1). I know that least squares might not be optimal in this context, but the models fit quite well. It is measured on a spatial data set and a model such as: sx9<-gnls(sprawlchngx9~ 1/(1+exp(-1*(beta0+ b1*hectare+ b2*discityge100kmean ))), data=data5, start=list(beta0 = 0,b1=0 , b2=0)) returns residuals with negative spatial autocorrelation at small lag distances between data points. The built in functions to handle spatially correlated residuals (corExp, etc) cannot (to my knowledge) handle negative spatial autocorrelation since their range parameters are bounded to be positive. To get around this issue, I extracted the residuals from the model above and computed Moran's I values at several lag distances. I then populated an NxN matrix with these values as determined by the lag distance between the row and column members. For instance, if observations 1 and 2 are 30 meters apart, and the Moran's I value for the lag distance between 0 and 50 m was computed as -0.50, then the matrix position (1,2) or (2,1) would be -0.50. I then passed this correlation matrix to gnls function with the following code: sx9<-gnls(sprawlchngx9~ 1/(1+exp(-1*(beta0+ b1*hectare+ b2*discityge100kmean ))), data=data5, start=list(beta0 = 0,b1=0 , b2=0, b16=0 , b17=0 , b18=0,b19=0,b20=0,b21=0 ), correlation=corSymm(corr9x[lower.tri(newcorr9x)],fixed=TRUE) ) where corr9x is the NxN correlation matrix that I created. Revisiting question 1, is there anything flawed in code or concept with this approach? Using residuals (), I have saved for this model "response", "pearson", "normalized" residuals. Based upon some digging, I thought that "normalized" would return residuals that are adjusted by my correlation matrix, but when variograms or correlograms are computed on all residuals they seem overly similar to those calculated on residuals NOT from the "normalized" option. Revisiting question 2, if residuals ( , type='normalized') the proper method to get residuals conditioned on the correlation structure (so that I can show whether or not the correlation= argument in gnls worked as intended)? If not, any suggestions? Thanks, Seth Myers [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] nearest neighbour using distance in metres
Dear Ross,
Euclidean Distance is not an units, it is a type of distance measurement an it
is always in input units!
The Euclidean Distance is the length of a rope between your hands that you hold
strained, but see here: http://en.wikipedia.org/wiki/Euclidean_distance
nndist(xySpatialTransformed@coords)
and maybe it is not very meaningful to calculate nndist with something else
than a planar coordinate system... or even equidistant coordinate system! (but
I'm not a cartographer)
Matteo
>>> Ross Ahmed 05/16/13 3:14 PM >>>
Hi all,
Following code computes Euclidian distance from each point to it nearest
neighbour:
library(spatstat)
y <- rnorm(30, 55.65, .008)
x <- rnorm(30, -1.82, .01)
xy <- cbind(x, y)
nndist(xy)
Is it possible to use distance in metres, not Euclidian distances? Here are
same points projected:
xySpatial <- SpatialPoints(xy, proj4string= CRS("+proj=longlat +datum=WGS84
+ellps=WGS84 +towgs84=0,0,0"))
xySpatialTransformed <- spTransform(xySpatial, CRS("+proj=laea +lon_0=-1.92
+lat_0=55.7 +datum=WGS84"))
Thanks
Ross
[[alternative HTML version deleted]]
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] nearest neighbour using distance in metres
Hi all,
Following code computes Euclidian distance from each point to it nearest
neighbour:
library(spatstat)
y <- rnorm(30, 55.65, .008)
x <- rnorm(30, -1.82, .01)
xy <- cbind(x, y)
nndist(xy)
Is it possible to use distance in metres, not Euclidian distances? Here are
same points projected:
xySpatial <- SpatialPoints(xy, proj4string= CRS("+proj=longlat +datum=WGS84
+ellps=WGS84 +towgs84=0,0,0"))
xySpatialTransformed <- spTransform(xySpatial, CRS("+proj=laea +lon_0=-1.92
+lat_0=55.7 +datum=WGS84"))
Thanks
Ross
[[alternative HTML version deleted]]
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
