[R-sig-Geo] R-stats - GRASS6

Roger Bivand Roger.Bivand at nhh.no
Wed Jun 1 13:33:26 CEST 2005


On Tue, 31 May 2005, Markus Neteler wrote:

> Dear Roger,
> 
> On Tue, May 31, 2005 at 06:53:31PM +0200, Roger Bivand wrote:
> > On Tue, 31 May 2005, Markus Neteler wrote:
> > 
> > > Dear list members,
> > > 
> > > I am currently trying out the new GRASS 6/R-stats interface.
> > > I installed all software today (R-2.1.0, all interfaces, using
> > > GRASS 6.1-CVS).
> > > 
> > > There is a little problem which I don't understand:
> > > 
> > > > library(spgrass6)
> > > > G <- gmeta6()
> > > > str(G)
> > > List of 28
> > >  $ GISDBASE         : chr "/ssi0/ssi/neteler/grassdata"
> > >  $ LOCATION_NAME    : chr "spearfish60"
> > >  $ MAPSET           : chr "neteler"
> > >  $ GRASS_DB_ENCODING: chr "utf-8"
> > >  $ DEBUG            : chr "0"
> > >  $ MONITOR          : chr "x0"
> > >  $ GRASS_GUI        : chr "text"
> > >  $ projection       : chr "1 (UTM)"
> > > [...]
> > > 
> > > > system("r.info -t landcover.orig")
> > > datatype=CELL
> > > > landcover <- readCELL6sp("landcover.orig", cat = TRUE)
> > > Error in "names<-.default"(`*tmp*`, value = "landcover.orig") :
> > >         'names' attribute [1] must be the same length as the vector [0]
> > 
> > I can't replicate this, I'm afraid (R 2.1.0, Linux, sp 0.7-7, spgrass6
> > 0.1-4 (from CVS on sourceforge, GRASS 6.0.0 release). I don't think the sp 
> > or spgrass6 versions make a difference. 
> > 
> > Could you say debug(readCELL6sp) and re-run this please - that should show 
> > which line is failing. Once we know that, please check the arguments to 
> > that function for sanity. I don't think it is important, but I'm running 
> > under LANG=en_GB ; export LANG, not UTF8. My current guess is that the 
> > problem is in:
> > 
> > res <- read.asciigrid(tmpfl, colname = vname, proj4string = p4)
> > 
> > where vname is the layer name. If you additionally do 
> > debug(read.asciigrid), it may fail at: names(df) = colname - please see 
> > what is in colname, and if length(colname) is the same as length(df). One 
> > possibility is that there is something amiss in read.asciigrid().
> > 
> > The landcover.orig layer is strange because there are no category labels 
> > at all, so:
> 
> [... took the wrong landuse map, but it doesn't change ...]
> 
> > library(spgrass6)
> Loading required package: rgdal
> Loading required package: abind
> Loading required package: pixmap
> Geospatial Data Abstraction Library extensions to R successfully loaded
> Loading required package: sp
> > G <- gmeta6()
> 
> > landcover <- readCELL6sp("landcover.30m", cat = TRUE)
> Error in "names<-.default"(`*tmp*`, value = "landcover.30m") :
>         'names' attribute [1] must be the same length as the vector [0]
> > debug(readCELL6sp)
> 
> > system("r.stats -l -q landcover.30m | head -20")
> 11 Open Water
> 21 Low Intensity Residential
> 22 High Intensity Residential
> 23 Commercial/Industrial/Transportation
> 31 Bare Rock/Sand/Clay
> 32 Quarries/Strip Mines/Gravel Pits
> 41 Deciduous Forest
> 42 Evergreen Forest
> 43 Mixed Forest
> 51 Shrubland
> 71 Grasslands/Herbaceous
> 81 Pasture/Hay
> 82 Row Crops
> 83 Small Grains
> 85 Urban/Recreational Grasses
> 91 Woody Wetlands
> 92 Emergent Herbaceous Wetlands
> * no data
> 
> > landcover <- readCELL6sp("landcover.30m", cat = TRUE)
> debugging in: readCELL6sp("landcover.30m", cat = TRUE)
> debug: {
> [...]
> }
> Browse[1]>
> debug: tmpfl <- tempfile()
> Browse[1]>
> debug: system(paste("r.out.arc input=", vname, " output=", tmpfl, sep = ""))
> Browse[1]>
> debug: library(sp)
> Browse[1]>
> debug: p4 <- CRS(system("g.proj -j -f", intern = TRUE))
> Browse[1]>
> debug: res <- read.asciigrid(tmpfl, colname = vname, proj4string = p4)
> Browse[1]>
> debug: if (cat) {
>     cats <- strsplit(system(paste("r.stats -l -q", vname), intern = TRUE),
>         " ")
>     catnos <- sapply(cats, function(x) x[1])
>     catlabs <- sapply(cats, function(x) paste(x[-1], collapse = " "))
>     if (any(!is.na(match(catnos, "*")))) {
>         isNA <- which(catnos == "*")
>         catnos <- catnos[-isNA]
>         catlabs <- catlabs[-isNA]
>     }
>     res at data[, 1] <- factor(res at data[, 1], levels = catnos, labels = catlabs)
> } else {
>     res at data[, 1] <- as.integer(res at data[, 1])
> }
> Browse[1]>
> debug: cats <- strsplit(system(paste("r.stats -l -q", vname), intern = TRUE),
>     " ")
> Browse[1]>
> debug: catnos <- sapply(cats, function(x) x[1])
> Browse[1]>
> debug: catlabs <- sapply(cats, function(x) paste(x[-1], collapse = " "))
> Browse[1]>
> debug: if (any(!is.na(match(catnos, "*")))) {
>     isNA <- which(catnos == "*")
>     catnos <- catnos[-isNA]
>     catlabs <- catlabs[-isNA]
> }
> Browse[1]>
> debug: isNA <- which(catnos == "*")
> Browse[1]>
> debug: catnos <- catnos[-isNA]
> Browse[1]>
> debug: catlabs <- catlabs[-isNA]
> Browse[1]>
> debug: res at data[, 1] <- factor(res at data[, 1], levels = catnos, labels = catlabs)
> Browse[1]>
> Error in "names<-.default"(`*tmp*`, value = "landcover.30m") :
>         'names' attribute [1] must be the same length as the vector [0]
> >  
> 
> 
> > debug(read.asciigrid)
> > undebug(readCELL6sp)
> 
> > landcover <- readCELL6sp("landcover.30m", cat = TRUE)
> [...]
> Browse[1]>
> debug: df = data.frame(map)
> Browse[1]>
> debug: names(df) = colname
> Browse[1]> length(colname)
> [1] 1
> Browse[1]> length(df)
> [1] 1
> Browse[1]> colname
> [1] "landcover.30m"
> Browse[1]> df
> [...all values...]
> 
> > str(df)
> `data.frame':   302418 obs. of  1 variable:
>  $ map: num  NA NA NA NA NA NA NA NA NA NA ...
> Browse[1]>
> debug: grid = GridTopology(c(xllcenter, yllcenter), rep(cellsize, 2),
>     c(ncols, nrows))
> Browse[1]>
> debug: SpatialGridDataFrame(grid, data = df, proj4string = proj4string)
> Browse[1]>
> exiting from: read.asciigrid(tmpfl, colname = vname, proj4string = p4)
> Error in "names<-.default"(`*tmp*`, value = "landcover.30m") :
>         'names' attribute [1] must be the same length as the vector [0]
> >    
> 
> I don't know if this is already helpful.

Yes, looks as though the changes we made to sp in Valencia (you have sp 
0.7-7, right?) are getting in the way. When the data is being made into a 
factor, the name of the data seems to be getting messed up. I have:

readCELL6sp <- function(vname, cat=FALSE) {
	tmpfl <- tempfile()
	system(paste("r.out.arc input=", vname, " output=", tmpfl, sep=""))
	library(sp)
	p4 <- CRS(system("g.proj -j -f", intern=TRUE))
	res <- read.asciigrid(tmpfl, colname=vname, proj4string=p4)
	if (cat) {
		cats <- strsplit(system(paste("r.stats -l -q", vname), 
			intern=TRUE), " ")
		catnos <- sapply(cats, function(x) x[1])
		catlabs <- sapply(cats, function(x) paste(x[-1], collapse=" "))
		if (any(!is.na(match(catnos, "*")))) {
			isNA <- which(catnos == "*")
			catnos <- catnos[-isNA]
			catlabs <- catlabs[-isNA]
		}
		res at data[[1]] <- factor(res at data[[1]], levels=catnos, labels=catlabs)
	} else {
		res at data[[1]] <- as.integer(res at data[[1]])
	}
	res
}

on my machine, I think it is on CVS, but cannot check now. The change 
Edzer and I made to sp was the same that was made to the R/GRASS interface 
years ago, that is do not use data frames where possible, because they 
generate unique character row names, which we do not need. So we changed 
from the @data slot being a data frame that can be accessed like a matrix, 
to a list of equal length vectors, which can only be accessed by list 
member. Please try this change in the function. The error message isn't 
the same in the following, but you can see what is going on (unless I'm 
still wrong):

> df <- data.frame(a=sample(1:5,100,replace=TRUE))
> str(df)
`data.frame':   100 obs. of  1 variable:
 $ a: int  3 2 2 2 3 1 3 1 5 2 ...
> ldf <- as.list(df)
> str(ldf)
List of 1
 $ a: int [1:100] 3 2 2 2 3 1 3 1 5 2 ...
> df1 <- df
> df1[,1] <- factor(df[,1])
> str(df1)
`data.frame':   100 obs. of  1 variable:
 $ a: Factor w/ 5 levels "1","2","3","4",..: 3 2 2 2 3 1 3 1 5 2 ...
> ldf1 <- ldf
> ldf1[,1] <- factor(ldf[,1])
Error in ldf[, 1] : incorrect number of dimensions
> ldf1[[1]] <- factor(ldf[[1]])
> str(ldf1)
List of 1
 $ a: Factor w/ 5 levels "1","2","3","4",..: 3 2 2 2 3 1 3 1 5 2 ...

Best wishes,

Roger

PS. I think 0.7-7 read.asciigrid() is also wrong, that is it should not 
use data.frame() anyway, but I can't check that now.

> 
> Best regards
> 
>  Markus
> 

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no




More information about the R-sig-Geo mailing list