在mgcv中從vis.gam更改等高線圖中的顏色

我正在使用vis.gam()來繪製模型的結果。代碼如下：在mgcv中從vis.gam更改等高線圖中的顏色

vis.gam(model,view=c("CHLA","SST"),xlim=c(0,15),ylim=c(10,30),xlab="CHL-a (mg.m-3)",ylab="SST (°C)",plot.type="contour",type="response",color="terrain",n.grid=200)

我想在下面的代碼中提到，從我自己的調色板設置新的顏色：

# colors palette definition 
    jet.colors <-colorRampPalette(c("#00007F", "blue", "#007FFF", "cyan", 
          "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000")) 

    vis.gam(model,view=c("CHLA","SST"),xlim=c(0,15),ylim=c(10,30),xlab="CHL-a (mg.m-3)",ylab="SST (°C)",plot.type="contour",type="response",color=jet.colors(10),n.grid=200)

不幸的是，我有這樣的錯誤：

vis.gam（model，view = c（「CHLA」，「SST」），xlim = c（0,15），ylim = c（10，： color scheme not recognized

我不知道如何使用jet.colors調色板功能vis.gam ...

感謝您的幫助！

來源

2014-02-13 habasque

從？vis.gam：'顏色：當se <= 0時用於圖的顏色方案。「topo」，「heat」，「cm」，「terrain」，「gray」或「bw」之一。方案「gray」和「bw」也修改se> 0時使用的顏色。所以顏色需要是一個字符串。 – EDi

您可以自定義vis.gam功能使用您的jet.colors斜坡：

下面是一個砍死版本（見我的在線評論）：

jet.colors <-colorRampPalette(c("#00007F", "blue", "#007FFF", "cyan", 
           "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000")) 

myvis.gam <- function (x, view = NULL, cond = list(), n.grid = 30, too.far = 0, 
      col = NA, color = "heat", contour.col = NULL, se = -1, type = "link", 
      plot.type = "persp", zlim = NULL, nCol = 50, ...) 
{ 
    fac.seq <- function(fac, n.grid) { 
    fn <- length(levels(fac)) 
    gn <- n.grid 
    if (fn > gn) 
     mf <- factor(levels(fac))[1:gn] 
    else { 
     ln <- floor(gn/fn) 
     mf <- rep(levels(fac)[fn], gn) 
     mf[1:(ln * fn)] <- rep(levels(fac), rep(ln, fn)) 
     mf <- factor(mf, levels = levels(fac)) 
    } 
    mf 
    } 
    dnm <- names(list(...)) 
    v.names <- names(x$var.summary) 
    if (is.null(view)) { 
    k <- 0 
    view <- rep("", 2) 
    for (i in 1:length(v.names)) { 
     ok <- TRUE 
     if (is.matrix(x$var.summary[[i]])) 
     ok <- FALSE 
     else if (is.factor(x$var.summary[[i]])) { 
     if (length(levels(x$var.summary[[i]])) <= 1) 
      ok <- FALSE 
     } 
     else { 
     if (length(unique(x$var.summary[[i]])) == 1) 
      ok <- FALSE 
     } 
     if (ok) { 
     k <- k + 1 
     view[k] <- v.names[i] 
     } 
     if (k == 2) 
     break 
    } 
    if (k < 2) 
     stop("Model does not seem to have enough terms to do anything useful") 
    } 
    else { 
    if (sum(view %in% v.names) != 2) 
     stop(paste(c("view variables must be one of", v.names), 
       collapse = ", ")) 
    for (i in 1:2) if (!inherits(x$var.summary[[view[i]]], 
           c("numeric", "factor"))) 
     stop("Don't know what to do with parametric terms that are not simple numeric or factor variables") 
    } 
    ok <- TRUE 
    for (i in 1:2) if (is.factor(x$var.summary[[view[i]]])) { 
    if (length(levels(x$var.summary[[view[i]]])) <= 1) 
     ok <- FALSE 
    } 
    else { 
    if (length(unique(x$var.summary[[view[i]]])) <= 1) 
     ok <- FALSE 
    } 
    if (!ok) 
    stop(paste("View variables must contain more than one value. view = c(", 
       view[1], ",", view[2], ").", sep = "")) 
    if (is.factor(x$var.summary[[view[1]]])) 
    m1 <- fac.seq(x$var.summary[[view[1]]], n.grid) 
    else { 
    r1 <- range(x$var.summary[[view[1]]]) 
    m1 <- seq(r1[1], r1[2], length = n.grid) 
    } 
    if (is.factor(x$var.summary[[view[2]]])) 
    m2 <- fac.seq(x$var.summary[[view[2]]], n.grid) 
    else { 
    r2 <- range(x$var.summary[[view[2]]]) 
    m2 <- seq(r2[1], r2[2], length = n.grid) 
    } 
    v1 <- rep(m1, n.grid) 
    v2 <- rep(m2, rep(n.grid, n.grid)) 
    newd <- data.frame(matrix(0, n.grid * n.grid, 0)) 
    for (i in 1:length(x$var.summary)) { 
    ma <- cond[[v.names[i]]] 
    if (is.null(ma)) { 
     ma <- x$var.summary[[i]] 
     if (is.numeric(ma)) 
     ma <- ma[2] 
    } 
    if (is.matrix(x$var.summary[[i]])) 
     newd[[i]] <- matrix(ma, n.grid * n.grid, ncol(x$var.summary[[i]]), 
          byrow = TRUE) 
    else newd[[i]] <- rep(ma, n.grid * n.grid) 
    } 
    names(newd) <- v.names 
    newd[[view[1]]] <- v1 
    newd[[view[2]]] <- v2 
    if (type == "link") 
    zlab <- paste("linear predictor") 
    else if (type == "response") 
    zlab <- type 
    else stop("type must be \"link\" or \"response\"") 
    fv <- predict.gam(x, newdata = newd, se.fit = TRUE, type = type) 
    z <- fv$fit 
    if (too.far > 0) { 
    ex.tf <- exclude.too.far(v1, v2, x$model[, view[1]], 
          x$model[, view[2]], dist = too.far) 
    fv$se.fit[ex.tf] <- fv$fit[ex.tf] <- NA 
    } 
    if (is.factor(m1)) { 
    m1 <- as.numeric(m1) 
    m1 <- seq(min(m1) - 0.5, max(m1) + 0.5, length = n.grid) 
    } 
    if (is.factor(m2)) { 
    m2 <- as.numeric(m2) 
    m2 <- seq(min(m1) - 0.5, max(m2) + 0.5, length = n.grid) 
    } 
    if (se <= 0) { 
    old.warn <- options(warn = -1) 
    av <- matrix(c(0.5, 0.5, rep(0, n.grid - 1)), n.grid, 
       n.grid - 1) 
    options(old.warn) 
    max.z <- max(z, na.rm = TRUE) 
    z[is.na(z)] <- max.z * 10000 
    z <- matrix(z, n.grid, n.grid) 
    surf.col <- t(av) %*% z %*% av 
    surf.col[surf.col > max.z * 2] <- NA 
    if (!is.null(zlim)) { 
     if (length(zlim) != 2 || zlim[1] >= zlim[2]) 
     stop("Something wrong with zlim") 
     min.z <- zlim[1] 
     max.z <- zlim[2] 
    } 
    else { 
     min.z <- min(fv$fit, na.rm = TRUE) 
     max.z <- max(fv$fit, na.rm = TRUE) 
    } 
    surf.col <- surf.col - min.z 
    surf.col <- surf.col/(max.z - min.z) 
    surf.col <- round(surf.col * nCol) 
    con.col <- 1 
    if (color == "heat") { 
     pal <- heat.colors(nCol) 
     con.col <- 3 
    } 
    else if (color == "topo") { 
     pal <- topo.colors(nCol) 
     con.col <- 2 
    } 
    else if (color == "cm") { 
     pal <- cm.colors(nCol) 
     con.col <- 1 
    } 
    else if (color == "terrain") { 
     pal <- terrain.colors(nCol) 
     con.col <- 2 
    } 
    else if (color == "gray" || color == "bw") { 
     pal <- gray(seq(0.1, 0.9, length = nCol)) 
     con.col <- 1 
    } 
    ### customized here 
    else if (color == 'jet') { 
     pal <- jet.colors(nCol) 
     con.col = 1 
    } 
    #### 
    else stop("color scheme not recognised") 
    if (is.null(contour.col)) 
     contour.col <- con.col 
    surf.col[surf.col < 1] <- 1 
    surf.col[surf.col > nCol] <- nCol 
    if (is.na(col)) 
     col <- pal[as.array(surf.col)] 
    z <- matrix(fv$fit, n.grid, n.grid) 
    if (plot.type == "contour") { 
     stub <- paste(ifelse("xlab" %in% dnm, "", ",xlab=view[1]"), 
        ifelse("ylab" %in% dnm, "", ",ylab=view[2]"), 
        ifelse("main" %in% dnm, "", ",main=zlab"), ",...)", 
        sep = "") 
     if (color != "bw") { 
     txt <- paste("image(m1,m2,z,col=pal,zlim=c(min.z,max.z)", 
        stub, sep = "") 
     eval(parse(text = txt)) 
     txt <- paste("contour(m1,m2,z,col=contour.col,zlim=c(min.z,max.z)", 
        ifelse("add" %in% dnm, "", ",add=TRUE"), ",...)", 
        sep = "") 
     eval(parse(text = txt)) 
     } 
     else { 
     txt <- paste("contour(m1,m2,z,col=1,zlim=c(min.z,max.z)", 
        stub, sep = "") 
     eval(parse(text = txt)) 
     } 
    } 
    else { 
     stub <- paste(ifelse("xlab" %in% dnm, "", ",xlab=view[1]"), 
        ifelse("ylab" %in% dnm, "", ",ylab=view[2]"), 
        ifelse("main" %in% dnm, "", ",zlab=zlab"), ",...)", 
        sep = "") 
     if (color == "bw") { 
     op <- par(bg = "white") 
     txt <- paste("persp(m1,m2,z,col=\"white\",zlim=c(min.z,max.z) ", 
        stub, sep = "") 
     eval(parse(text = txt)) 
     par(op) 
     } 
     else { 
     txt <- paste("persp(m1,m2,z,col=col,zlim=c(min.z,max.z)", 
        stub, sep = "") 
     eval(parse(text = txt)) 
     } 
    } 
    } 
    else { 
    if (color == "bw" || color == "gray") { 
     subs <- paste("grey are +/-", se, "s.e.") 
     lo.col <- "gray" 
     hi.col <- "gray" 
    } 
    else { 
     subs <- paste("red/green are +/-", se, "s.e.") 
     lo.col <- "green" 
     hi.col <- "red" 
    } 
    if (!is.null(zlim)) { 
     if (length(zlim) != 2 || zlim[1] >= zlim[2]) 
     stop("Something wrong with zlim") 
     min.z <- zlim[1] 
     max.z <- zlim[2] 
    } 
    else { 
     z.max <- max(fv$fit + fv$se.fit * se, na.rm = TRUE) 
     z.min <- min(fv$fit - fv$se.fit * se, na.rm = TRUE) 
    } 
    zlim <- c(z.min, z.max) 
    z <- fv$fit - fv$se.fit * se 
    z <- matrix(z, n.grid, n.grid) 
    if (plot.type == "contour") 
     warning("sorry no option for contouring with errors: try plot.gam") 
    stub <- paste(ifelse("xlab" %in% dnm, "", ",xlab=view[1]"), 
        ifelse("ylab" %in% dnm, "", ",ylab=view[2]"), ifelse("zlab" %in% 
                     dnm, "", ",zlab=zlab"), ifelse("sub" %in% dnm, 
                             "", ",sub=subs"), ",...)", sep = "") 
    txt <- paste("persp(m1,m2,z,col=col,zlim=zlim", ifelse("border" %in% 
                  dnm, "", ",border=lo.col"), stub, sep = "") 
    eval(parse(text = txt)) 
    par(new = TRUE) 
    z <- fv$fit 
    z <- matrix(z, n.grid, n.grid) 
    txt <- paste("persp(m1,m2,z,col=col,zlim=zlim", ifelse("border" %in% 
                  dnm, "", ",border=\"black\""), stub, sep = "") 
    eval(parse(text = txt)) 
    par(new = TRUE) 
    z <- fv$fit + se * fv$se.fit 
    z <- matrix(z, n.grid, n.grid) 
    txt <- paste("persp(m1,m2,z,col=col,zlim=zlim", ifelse("border" %in% 
                  dnm, "", ",border=hi.col"), stub, sep = "") 
    eval(parse(text = txt)) 
    } 
}

對於使用只需使用color = 'jet'：

library(mgcv) 
set.seed(0) 
n<-200;sig2<-4 
x0 <- runif(n, 0, 1);x1 <- runif(n, 0, 1) 
x2 <- runif(n, 0, 1) 
y<-x0^2+x1*x2 +runif(n,-0.3,0.3) 
g<-gam(y~s(x0,x1,x2)) 
myvis.gam(g,ticktype="detailed",color='heat',theta=-35) 
myvis.gam(g,ticktype="detailed",color='jet',theta=-35)

enter image description here

來源

2014-02-13 10:01:28 EDi

在mgcv中從vis.gam更改等高線圖中的顏色

回答

相關問題