+ Show code- Hide code
library(OpasnetUtils)
ograph <- function( # Määritellään yleisfunktio peruskuvaajan piirtämiseen.
ovariable,
x,
y = character(),
type = character(),
other = character(),
fill = NA,
...
) {
if(class(ovariable) == "ovariable") {
if(nrow(ovariable@output) == 0) ovariable <- EvalOutput(ovariable)
data <- ovariable@output
title <- ovariable@name
if(length(y) == 0) y <- paste(title, "Result", sep = "")
} else {
data <- ovariable
title <- character()
if(length(y) == 0) y <- "Result"
}
if(length(type) == 0) {
if("Iter" %in% colnames(data)) type <- geom_boxplot() else type <- geom_bar(stat = "identity")
}
out <- ggplot(data, aes_string(x = x, y = y, fill = fill)) # Määritellään kuvan sarakkeet
out <- out + type
out <- out + theme_grey(base_size=24) # Fontin kokoa suurennetaan
out <- out + labs(
title = title,
y = paste(unique(data[[paste(title, "Yksikkö", sep = "")]]), sep = "", collapse = ", ")
)
out <- out + theme(axis.text.x = element_text(angle = 90, hjust = 1)) # X-akselin tekstit käännetään niin että mahtuvat
if(length(other) != 0) out <- out + other
return(out)
}
# fillna takes a data.frame and fills the cells with NA with each level in that column.
# object is the data.frame, marginals is a vector of columns (either column names or positions) that are to be filled.
# This version of fillna accepts column positions (as the previous version) and also column names in marginals.
fillna <- function (object, marginals) {
a <- dropall(object)
if(!is.numeric(marginals)) marginals <- match(marginals, colnames(object))
for (i in marginals) {
a[[i]] <- as.factor(a[[i]])
a1 <- a[!is.na(a[[i]]), ]
a2 <- a[is.na(a[[i]]), ][-i]
addition <- data.frame(A = levels(a[[i]]))
colnames(addition) <- colnames(a)[i]
a2 <- merge(addition, a2)
a <- rbind(a1, a2)
}
return(a)
}
## collapsemarg is a placeholder for a better functionality within CollapseMarginals.
## It takes an ovariable, and summarises all indices in cols using tapply and a user-defined function.
## However, you can also use function "pick" to select locations defined in a list picks found in indices cols.
## Function "unkeep" simply drops the unkept indices without any other operation.
## The output is an ovariable with the same name as the input.
## This was first created for [[:op_fi:Radonin terveysvaikutukset]]
collapsemarg <- function(variable, cols, fun = "sum", picks = list(), ...) { # cols is a character vector, while probs is a list
out <- dropall(variable@output)
marginals <- colnames(out)[variable@marginal]
if(tolower(fun) == "unkeep") { # The function must be a string, otherwise this row will fail.
out <- out[!colnames(out) %in% cols]
} else {
if(tolower(fun) == "pick") {
for(i in cols) {
out <- out[out[[i]] %in% picks[[match(i, cols)]] , ]
}
cols <- "" # Those locations that were picked are still marginals.
} else {
margtemp <- colnames(out)[colnames(out) %in% marginals & !colnames(out) %in% cols]
# You must leave at least one index.
out <- as.data.frame(as.table(tapply(result(variable), out[margtemp], fun)))
out <- out[!is.na(out$Freq) , ]
colnames(out)[colnames(out) == "Freq"] <- ifelse(
length(variable@name) == 0,
"Result",
paste(variable@name, "Result", sep = "")
)
}
}
variable@output <- out
variable@marginal <- colnames(out) %in% marginals & ! colnames(out) %in% cols
return(variable)
}
MyPointKML <- function( # The function creates a KML fille from a SpatialPointsDataFrame object.
obj = NULL, # Spatial object with the data. A SpatialPointsDataFrame.
kmlname = "", # Name of the KML fille (does this show on the map?)
kmldescription = "", # Description of the KML fille (does this show on the map?)
name = NULL, # Name for each datapoint (vector with the same length as data in obj).
description = "", # Descrtion of each datapoint (vector with the same length as data in obj).
icon = "http://maps.google.com/mapfiles/kml/pal4/icon24.png", # Icon shown on pin (?)
col=NULL # I don't know what this does.
) {
cat("This function MyPointKML is depreciated. Use google.point_kml in OpasnetUtilsExt instead.\n")
if (is.null(obj))
return(list(header = c("<?xml version=\"1.0\" encoding=\"UTF-8\"?>",
"<kml xmlns=\"http://earth.google.com/kml/2.2\">",
"<Document>", paste("<name>", kmlname, "</name>",
sep = ""), paste("<description><![CDATA[", kmldescription,
"]]></description>", sep = "")), footer = c("</Document>",
"</kml>")))
if (class(obj) != "SpatialPointsDataFrame")
stop("obj must be of class 'SpatialPointsDataFrame' [package 'sp']")
if (is.null(name)) {
name = c()
for (i in 1:nrow(obj)) name <- append(name, paste("site",
i))
}
if (length(name) < nrow(obj)) {
if (length(name) > 1)
warning("kmlPoints: length(name) does not match nrow(obj). The first name will be replicated.")
name <- rep(name, nrow(obj))
}
if (length(description) < nrow(obj)) {
if (length(description) > 1)
warning("kmlPoints: length(description) does not match nrow(obj). The first description will be replicated.")
description <- rep(description, nrow(obj))
}
if (length(icon) < nrow(obj)) {
if (length(icon) > 1)
warning("kmlPoints: length(icon) does not match nrow(obj). Only the first one will be used.")
icon <- icon[1]
}
# This is some kind of a colour definition
col2kmlcolor <- function(col)
paste(rev(sapply(
col2rgb(col, TRUE),
function(x) sprintf("%02x", x))
), collapse = "")
kml <- kmlStyle <- ""
# Create the KML fille.
kmlHeader <- c("<?xml version=\"1.0\" encoding=\"UTF-8\"?>","<kml xmlns=\"http://earth.google.com/kml/2.2\">", "<Document>")
kmlFooter <- c("</Document>", "</kml>")
# Create rows to the KML fille from data in obj.
for (i in 1:nrow(obj)) {
point <- obj[i, ]
pt_style <- paste("#style", ifelse(length(icon) == 1, 1, i), sep = "")
kml <- append(kml, "<Placemark>")
kml <- append(kml, paste(
" <description><![CDATA[",
name[i],
": ",
description[i],
"]]></description>",
sep = ""
))
#kml <- append(kml, "<Style><IconStyle>")
#kml <- append(kml, paste("<color>", col2kmlcolor(col[i]), "</color>", sep =""))
#kml <- append(kml, paste(" <Icon><href>", icon, "</href></Icon>", sep = ""))
#kml <- append(kml, "<scale>0.300000</scale>")
#kml <- append(kml, "</IconStyle></Style>")
kml <- append(kml, " <Point>")
kml <- append(kml, " <coordinates>")
kml <- append(kml, paste(point@coords[1], point@coords[2], sep = ","))
kml <- append(kml, " </coordinates>")
kml <- append(kml, " </Point>")
kml <- append(kml, "</Placemark>")
}
return(paste(paste(c(kmlHeader, kmlStyle, kml, kmlFooter), sep = "", collapse = "\n"), collapse="\n", sep = ""))
}
ova2spat <- function( # This function converts an ovariable into a SpatialPointsDataFrame.
ovariable, # An evaluated ovariable that has coordinate indices.
coords, # The names of the coordinate indices as a character vector, first x then y.
proj4string # Projection identifier or specification as character string. See http://spatialreference.org/
) {
temp <- ovariable@output
# Transform coordinates into numeric format.
for(i in coords) {
if(is(temp[[i]], "factor")) temp[[i]] <- levels(temp[[i]])[temp[[i]]]
if(is(temp[[i]], "character")) temp[[i]] <- as.numeric(temp[[i]])
}
# Define the coordinate points first, then add other ovariable output to it.
sp <- SpatialPoints(temp[coords], CRS(proj4string))
out <- SpatialPointsDataFrame(sp, temp[!colnames(temp) %in% coords])
#Transform the projection to longitude-latitude system.
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
out <- spTransform(out,epsg4326String)
return(out)
}
# MyRmap is a function for creating static Google maps as png.
# It is based on MyMap function without the "file destination" parameter
# Requires RgoogleMaps package
MyRmap <- function (
shp, # a spatial data object
plotvar, # Name of the column that has the values to be illustrated on the map
pch = 19, # Shape of the point (19: circle)
cex = 0.3, # Size of the point
legend_title = "", # Title of the legend
legend_position = "topleft",
numbins = 8, # Number of colour bins in graph
center, # center of the map
size = c(640, 480), # size of the map. This produces the right dimensions in Opasnet.
MINIMUMSIZE = FALSE,
RETURNIMAGE = TRUE,
GRAYSCALE = FALSE,
NEWMAP = TRUE,
zoom,
verbose = 1,
...
) {
plotvar <- shp[[plotvar]]
plotclr <- brewer.pal(numbins, "Spectral")
classes <- classIntervals(plotvar, numbins, style = "quantile")
colcode <- findColours(classes, plotclr)
latR <- shp@coords[ , 2]
lonR <- shp@coords[ , 1]
#get the bounding box:
bb <- qbbox(lat = latR, lon = lonR)
if (missing(zoom))
zoom <- min(MaxZoom(latR, lonR, size))
if (missing(center)) {
lat.center <- mean(latR)
lon.center <- mean(lonR)
}
else {
lat.center <- center[1]
lon.center <- center[2]
}
if (MINIMUMSIZE) {
ll <- LatLon2XY(latR[1], lonR[1], zoom) # I think the latR and lonR are used here differently than how they
ur <- LatLon2XY(latR[2], lonR[2], zoom) # are used elsewhere. Thus, if MINIMUMSIZE = TRUE, you may see problems.
cr <- LatLon2XY(lat.center, lon.center, zoom)
ll.Rcoords <- Tile2R(ll, cr)
ur.Rcoords <- Tile2R(ur, cr)
if (verbose > 1) {
cat("ll:")
print(ll)
print(ll.Rcoords)
cat("ur:")
print(ur)
print(ur.Rcoords)
cat("cr:")
print(cr)
}
size[1] <- 2 * max(c(ceiling(abs(ll.Rcoords$X)), ceiling(abs(ur.Rcoords$X)))) + 1
size[2] <- 2 * max(c(ceiling(abs(ll.Rcoords$Y)), ceiling(abs(ur.Rcoords$Y)))) + 1
if (verbose) cat("new size: ", size, "\n")
}
MyMap <- GetMap(
center = c(lat.center, lon.center),
zoom = zoom,
size = size,
RETURNIMAGE = RETURNIMAGE,
GRAYSCALE = GRAYSCALE,
verbose = verbose,
...
)
PlotOnStaticMap(MyMap) # Plot an empty map.
PlotOnStaticMap( # Plot the data points on the map.
MyMap,
lat = latR,
lon = lonR,
pch = pch,
cex = cex,
col = colcode,
add = T
)
legend( # Plot the legend on the map.
legend_position,
legend = names(attr(colcode, "table")),
title = legend_title,
fill = attr(colcode, "palette"),
cex = 1.0,
bty = "y",
bg = "white"
)
}
MyPlotKML <- function(
shp, # a SpatialPointDataFrame object.
result = "Result", # The name of result column in shp.
rasterization = TRUE, # Whether to rasterize the data or not.
ncols = 32, # Number or columns in the raster.
nrows = 32, # Number of rows in the raster.
fun = mean # function to aggregate data points to the raster.
) {
cat("Consider merging this function MyPolotKML with google.show_raster_on_maps in OpasnetUtilsExt.\n")
if(rasterization) {
#Create blank raster
rast <- raster()
#Set raster extent to that of point data
extent(rast) <-extent(shp)
#Choose number of columns and rows
ncol(rast) <- ncols
nrow(rast) <- nrows
#Rasterize point data
rast2 <- rasterize(shp, rast, shp[[result]], fun = fun)
}
start <- 0 # min(shp[[result]])
end <- max(shp[[result]])
steps <- approx(c(start,end),n=6)$y
colors <- rev(rainbow(length(steps), start=0, end=0.50))
# Create the colorstrip below the map.
par(mfrow=c(6,1), mar=c(3,1,0,1), cex = 1.5)
colorstrip <- function(colors, labels)
{
count <- length(colors)
image(
matrix(1:count, count, 1),
col = colors,
ylab = "",
axes = FALSE
)
axis(1,approx(c(0, 1), n=length(labels))$y, labels)
}
colorstrip(colors, steps)
#Plot data
google.show_raster_on_maps(rast2, col = colors, style = "height:500px;")
}
###### ovaMatrixProduct is a matrix multiplication (%*%) for ovariables.
ovaMatrixProduct <- function(
e1, # The first ovariable in the multiplication.
e2, # The second ovariable in the multiplication.
matrind.a, # The index names used in the matrix (row, column) for the first ovariable. Character vector, length 2.
matrind.b # The index names used in the matrix (row, column) for the second ovariable.
) {
#matrind.a <- c("A", "B")
#matrind.b <- c("B", "C")
#Take note of names and marginals
#Merge the two ovariables
#Make a data frame with all non-matrix indices, unique rows.
#Make a for loop and go through non-matrix index rows one by one.
## Slice the merged ovariable and take the result columns and the matrix indices for each original ovariable.
## Reshape the two new data.frames to wide format and convert to matrices.
## Make the matrix product.
## Convert back to data frame.
## Append to previous result data.frame together with non-matrix indices.
# Return a new ovariable.
matricise <- function( # Produce a matrix of the key data from a data.frame.
x, # data.frame with the data for making a matrix.
matrind, # Index names in x for the matrix.
name # Name of the result column in x.
) {
out <- unique(x[c(matrind, name)])
out <- reshape(out, idvar = matrind[1], timevar = matrind[2], direction = "wide")
colnames(out) <- gsub(paste(name, ".", sep = ""), "", colnames(out))
rownames(out) <- out[[1]]
out[[1]] <- NULL
out <- as.matrix(out)
return(out)
}
temp <- merge(e1, e2) # Combine ovariables to make all indices match in the future calculations.
# Create a data.frame that contains unique rows consisting of other marginal indices than those needed in the matrices.
temp2 <- c(colnames(e1@output)[e1@marginal], colnames(e2@output)[e2@marginal])
temp2 <- temp2[! temp2 %in% c(matrind.a, matrind.b)]
nonmatrixind <- unique(temp@output[temp2])
# Create an output data.frame that collects a matrix product result for each rows of nonmatrixind.
out <- data.frame()
for(i in 1:nrow(nonmatrixind))
{
temp3 <- merge(temp@output, nonmatrixind[i, , drop = FALSE])
# Create matrices for ovariables a and b from data.frame temp3.
tempa <- matricise(temp3, matrind.a, comment(result(e1)))
tempb <- matricise(temp3, matrind.b, comment(result(e2)))
matr <- tempa %*% tempb # The actual matrix product.
# Convert matrix back to data.frame and combine with other results.
matr <- as.data.frame(as.table(matr))
colnames(matr) <- c(matrind.a[1], matrind.b[2], "Result")
matr <- cbind(nonmatrixind[i, , drop = FALSE], matr)
out <- rbind(out, matr)
}
out <- new("ovariable", output = out)
out <- CheckMarginals(out, deps = list(e1, e2), verbose = FALSE)
return(out)
}
# Merge all but show_bins largest bins of indices cols to 'Other'.
truncateIndex <- function( # Truncates an index to contain only the largest index bins.
obj, # ovariable to use.
cols, # names of the columns to truncate.
bins = 10# Number of bins to show. Other bins will be merged to "Other".
) {
for(i in 1:length(cols))
{
sums <- as.data.frame(as.table(tapply(result(obj), obj@output[cols[i]], sum)))
bins <- min(bins, length(levels(obj@output[[cols[i]]])))
limit <- sort(sums$Freq, decreasing = TRUE)[bins]
keeps <- sums[sums$Freq > limit , cols[i]]
keeps <- levels(keeps)[keeps]
levels(obj@output[[cols[i]]]) <- ifelse(
levels(obj@output[[cols[i]]]) %in% keeps,
levels(obj@output[[cols[i]]]),
"Other"
)
}
return(obj)
}
objects.store(ograph, fillna, collapsemarg, MyPointKML, ova2spat, MyRmap, MyPlotKML, ovaMatrixProduct, truncateIndex)
cat("The following objects are stored: ograph, fillna, collapsemarg, MyPointKML, ova2spat, MyRmap, MyPlotKML, ovaMatrixProduct, truncateIndex.\n")
| |
