[R] Unexpected behavior in PBSmapping package

[D. Dailey]

Thank you, Rowan, for the explanation and another workaround.

I think my major confusion stemmed from the idea that the plotted point 
colors would change only as a result of changing the limits of the plot. 
  I was and still am confused by the disparity between the plot using 
PBSmapping and the following, using standard graphics, where the colors 
of the points do not depend on the limits of the plot:

## Begin Example ##

X <- c( 494,497,500,503)
Y <- 1000 + X
colors <- c( 'red', 'blue', 'green', 'yellow' )

par( mfrow = c(2,1) )
plot( X, Y, pch=16, col=colors, cex = 3 )
text( X, Y, labels = toupper( substr( colors, 1, 1 ) ), font=2 )

plot( X, Y, pch=16, col=colors, cex = 3,
   xlim=c( 499, 505 ), ylim=c( 1499, 1505 ) )
text( X, Y, labels = toupper( substr( colors, 1, 1 ) ), font=2 )

## End Example ##

If I understand PBSmapping correctly, I think I see that plotPoints 
first takes a subset of the data argument, keeping only those points 
that lie within the figure's limits, and then applies the vectors of 
parameters such as col (probably also cex and pch though I haven't 
checked them), even though the col vector almost certainly now is 
mismatched with the data.  It sounds like you're saying that if these 
parameters are defined by polyProps instead, a similar subset of 
polyProps is taken before accessing the col (and cex, pch) columns.

In my toy example, the color was a column of the data frame.  In my 
actual work, I've been handling the color assignment outside the data 
frame and sending the color vector as an argument to plotPoints 
(actually, to addPoints, but the problem is the same).

So if I want to use col as an argument directly, it seems that I have to 
know which points are going to fall outside the limits of my plot.  This 
seems an unreasonable burden to place on the user.

In the documentation for plotPoints(), col is defined as a "vector of 
colours (cycled by EID or PID)"... but the behavior is more as though it 
were a "vector of colours (cycled by EID or PID *among the subset of 
points remaining in the plotting region*)".  Your reply seems to suggest 
that this latter was the intended behavior.

My biggest concern is that, using the col argument, the colors assigned 
to the points are dependent on the limits of the plotting region.  As a 
further thought, consider the following adapted from the example code 
for plotPoints:

## Begin example ##

require( PBSmapping )
data( nepacLL )
data( surveyData )

par( mfrow=c(2,1) )
# Example code from plotPoints:
plotMap(nepacLL, xlim=c(-136, -125), ylim=c(48, 57))
addPoints(surveyData, col=1:7, pch=16)

# Now "zoom in" to latitudes 52 to 57 degrees north
plotMap(nepacLL, xlim=c(-136, -125), ylim=c(52, 57))
addPoints(surveyData, col=1:7, pch=16)

## End example ##

The points in Dixon Entrance (along the north edge of Graham Island) 
have different colors in the second plot than they do in the first.  And 
they could have even different colors if the plot had different xlim and 
ylim values.  I hope you can help me understand why this is desirable 
behavior of these functions.


I do have to say that I greatly appreciate the work and thought that has 
gone into, and continues to go into, the development of this package. 
The package has enabled me to provide a great service to a local school 
district that was struggling with the closure of a couple of schools and 
the subsequent redistribution of the students among the remaining buildings.

Again, thanks for the ideas for the workaround.  I hope my comments will 
provide good substance for further discussions within your workgroup.

--David Dailey
Shoreline, Washington, USA

On 8/13/2007 10:24 AM, Haigh, Rowan wrote:
> Hi David,
> Thanks for the message. Essentially, the argument 'col' is treated like
> a vector which remains independent from (X,Y) and is simply repeated as
> need for the number of points. If you specify
> col=c("red","blue","green","yellow") as below, and there are only 2
> points, only "red" and "blue" will be used (hence the perceived
> mismatch).
> 
> To keep your colour column matched with X,Y, there is an argument called
> polyProps that automatically looks for the column "col". Unfortunately,
> polyProps doesn't appear to be functional with EventData (bug, I think)
> but will work with PolyData. Take a look at the R-code below. I've
> redefined your EventData to be PolyData which can serve double duty as
> EventData and PolyData, i.e., you can pass "events" to plotPoints() as
> well as using "polyProps=events". I've also added the "label" field to
> events to simplify things, and used PBSmapping's addLabels() function. 
> 
> The plotPoints() code should have been able to do this with EventData
> also. We will have to revisit this code. Thanks David for taking the
> time to show us user-related problems. It helps us to re-evaluate
> functionalities. Undoubtedly, the package has many little imperfections
> but we try to iron them out over time.
> Cheers, Rowan
> 
> ### Begin Example ###
> library( PBSmapping )
> 
> # Define some EventData
> events <- as.PolyData( read.table( textConnection(
> 'PID X  Y  col
> 1 494 1494 red
> 2 497 1497 blue
> 3 500 1500 green
> 4 503 1503 yellow' ), header=TRUE, strings=FALSE ),
> projection='UTM', zone=10 )
> events$label <- toupper(substring(events$col,1,1));
> 
> par( mfrow=c(3,1) )
> 
> # Plot the events with plot limits large enough to show
> # the full extent of all the symbols
> plotPoints( events, pch=16, cex=5, polyProps=events,
>    xlim=c(490,508), ylim=c(1490,1508),proj="UTM")
> addLabels(events, polyProps=events, font=2, cex=2, col=1) 
> 
> # Normal plot extents; partial symbols cut off by edges
> # of plotting region (as expected)
> plotPoints( events, pch=16, cex=5, polyProps=events,proj="UTM")
> addLabels(events, polyProps=events, font=2, cex=2, col=1) 
> 
> ## Now use more-restrictive plot limits
> plotPoints( events, pch=16, cex=5, polyProps=events,
>    xlim=c(499,505), ylim=c(1499,1505),proj="UTM")
> addLabels(events, polyProps=events, font=2, cex=2, col=1) 
> 
> ### End example ###
> 
> P.S. Once you've defined your EventData/PolyData/PolySet with attribute
> projection="UTM", plotMap() will attempt to use this but plotLines() and
> plotPoints() will simply tell you if it matched the specified projection
> in the function call. When you specified plotPoints(...,proj=TRUE), you
> are probably overriding the earlier specified UTM projection (of events)
> with proj=1. In this case, the 2 projections (UTM, 1) are the same. If
> your Poly dataset was set to projection="LL", they would not be.
> ======================================================
> 
> -----Original Message-----
> From: Schnute, Jon 
> Sent: Monday, August 13, 2007 8:19 AM
> To: Haigh, Rowan
> Subject: FW: Unexpected behavior in PBSmapping package
> 
> 
> Hi Rowan - Could you please check into this when you get a chance?
> Thanks - Jon 
> 
> -----Original Message-----
> From: D. Dailey [mailto:lists at CompassRoseEnterprises.com] 
> Sent: Thursday, August 09, 2007 11:38 PM
> To: r-help at stat.math.ethz.ch; Schnute, Jon
> Subject: Unexpected behavior in PBSmapping package
> 
> Using R 2.5.1 on Windows XP Professional, and PBSmapping package version
> 
> 2.51, I have encountered some behavior which puzzles me.  I am including
> 
> the package's listed maintainer on this email but also seek the thoughts
> 
> of the R-help community.
> 
> I have a set of EventData, which I want to plot as points, and to color 
> the points according to some criterion.  It turns out that some of my 
> points fall outside my desired plotting region.  It looks like this 
> causes the PBSmapping functions plotPoints and addPoints to incorrectly 
> deal with the color assignments.
> 
> Consider the following toy example:
> 
> ### Begin Example ###
> 
> library( PBSmapping )
> 
> # Define some EventData
> events <- as.EventData( read.table( textConnection(
> 'EID X  Y   Color
> 1 494 1494 red
> 2 497 1497 blue
> 3 500 1500 green
> 4 503 1503 yellow' ), header=TRUE, strings=FALSE ),
> proj='UTM', zone=10 )
> 
> par( mfrow=c(3,1) )
> 
> # Plot the events with plot limits large enough to show
> # the full extent of all the symbols
> plotPoints( events, pch=16, cex=5, col=events$Color,
>    xlim=c(490,508), ylim=c(1490,1508), proj=TRUE )
> with( events, text( X, Y, toupper( substr( Color, 1, 1 ) ),
>    font=2, cex=2 ) )
> 
> # Normal plot extents; partial symbols cut off by edges
> # of plotting region (as expected)
> plotPoints( events, pch=16, cex=5, col=events$Color, proj=TRUE )
> with( events, text( X, Y, toupper( substr( Color, 1, 1 ) ),
>    font=2, cex=2 ) )
> 
> ## Now use more-restrictive plot limits
> plotPoints( events, pch=16, cex=5, col=events$Color,
>    xlim=c(499,505), ylim=c(1499,1505), proj=TRUE )
> with( events, text( X, Y, toupper( substr( Color, 1, 1 ) ),
>    font=2, cex=2 ) )
> # Note that symbols are plotted in the right places (note text labels)
> # but colors are not as expected
> 
> ### End example ###
> 
> 
> For the moment, I have worked around this issue by using a "with( 
> events, points( ... ) )" construction, but this seems suboptimal; I 
> would prefer to use addPoints (which exhibits the same problem as 
> plotPoints does in the toy example above).  I would appreciate any 
> insights those on the list might have.
> 
> Please include me directly on any reply to the list, as I am at least a 
> couple weeks behind on reading the digested version of the list.  I see 
> that there have been no mentions of the PBSmapping package even in the 
> digests I have not yet read.
> 
> Session info:
> 
>  > sessionInfo()
> R version 2.5.1 (2007-06-27)
> i386-pc-mingw32
> 
> locale:
> LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
> States.1252;LC_MONETARY=English_United 
> States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
> 
> attached base packages:
> [1] "stats"     "graphics"  "grDevices" "utils"     "datasets" 
> "methods"   "base"
> 
> other attached packages:
> PBSmapping
>      "2.51"
> 
> 
> 
> --David Dailey
> Shoreline, Washington, USA