RSS Feed : rgraphgallery.blogspot.com

RSS Feed from http://rgraphgallery.blogspot.com

  • RG #110: 3D scatter plot with multiple series in Y axis
    X = seq(1, 100, 5)
    Y = seq (1, 100, 5)
    Z = rnorm (length (X), 10, 2)
    data1 <- data.frame (X, Y, )
    data2 <- data.frame (X, Y, Z1 = Z - 5)
    data3 <- data.frame (X, Y, Z1 = Z - 3)


    require(scatterplot3d)
    s3d <- scatterplot3d(data1, color = "blue", pch = 19, xlim=NULL, ylim=NULL, zlim= c(0, 20))
    s3d$points3d(data2, col = "red", pch = 18)
    s3d$points3d(data3, col = "green4", pch = 17)



     

    Read more »
  • RG#104: Rootogram plot
    require(latticeExtra)
    require(lattice)


    library(lattice)
    set.seed(1234)
    x1 <- rpois(1000, lambda = 25)


    plt <- rootogram(~x1, dfun = function(x) dpois(x, lambda = 25))
    plt

     
     

    Read more »
  • RG#109:small plot(s) with in a big plot
    require(ggplot2)
    library(gridBase)

    plot(cos, -pi, 2*pi, ylim = c(-1.3, 1.5), col = "red")
    myd <- data.frame (X = 1:10, Y = c(3, 4, 8, 7, 2, 1, 9, 4, 2, 3))
    qp <- qplot(X, Y, data=myd) + theme_bw()

    print(qp, vp=viewport(.65, .65, .25, .25))

     library(lattice)
    library(gridBase)

    plot.new()
    pushViewport(viewport())
    set.seed(1234)
    xvars <- rnorm(25, 5, 1)
    yvars <- rnorm(25, 5, 1)
    xyplot(yvars~xvars,  xlim = c(0, 10), ylim = c(0, 10) )
    pushViewport(viewport(x=.6,y=.85,width=.20,height=.15,just=c("left","top")))
    grid.rect()
    par(plt = gridPLT(), new=TRUE)
    plot(xvars,yvars)
    popViewport(2)






     
    Read more »
  • RG#107: Plot 3d horizontal lines (bars) over map (world and US example)
    library("maps")
    require(ggplot2)
    library(ggsubplot)

    world.map <- map("world", plot = FALSE, fill = TRUE)
    world_map <- map_data("world")
    require(lattice)
    require(latticeExtra)


     # Calculate the mean longitude and latitude per region (places where subplots are plotted)

    library(plyr)
    cntr <- ddply(world_map,.(region),summarize,long=mean(long),lat=mean(lat))



    # example data
     myd <- data.frame (region = c("USA","China","USSR","Brazil", "Australia","India", "Nepal", "Canada",
                                    "South Africa", "South Korea", "Philippines", "Mexico", "Finland",
                                     "Egypt", "Chile", "Greenland"),
                   frequency = c(501, 350, 233, 40, 350, 150, 180, 430, 233, 120, 96, 87, 340, 83, 99, 89))




    subsetcntr  <- subset(cntr, region %in% c("USA","China","USSR","Brazil", "Australia","India", "Nepal", "Canada",
                                    "South Africa", "South Korea", "Philippines", "Mexico", "Finland",
                                     "Egypt", "Chile", "Greenland"))


    simdat <- merge(subsetcntr, myd)
    colnames(simdat) <- c( "region","long","lat", "myvar" )



    panel.3dmap <- function(..., rot.mat, distance, xlim,
         ylim, zlim, xlim.scaled, ylim.scaled, zlim.scaled) {
           scaled.val <- function(x, original, scaled) {
          scaled[1] + (x - original[1]) * diff(scaled)/diff(original)
         }
           m <- ltransform3dto3d(rbind(scaled.val(world.map$x,
               xlim, xlim.scaled), scaled.val(world.map$y, ylim,
              ylim.scaled), zlim.scaled[1]), rot.mat, distance)
            panel.lines(m[1, ], m[2, ], col = "green4")
          }



    p2 <- cloud(myvar ~ long + lat, simdat, panel.3d.cloud = function(...) {
             panel.3dmap(...)
              panel.3dscatter(...)
     }, type = "h", col = "red", scales = list(draw = FALSE), zoom = 1.1,
                xlim = world.map$range[1:2], ylim = world.map$range[3:4],
              xlab = NULL, ylab = NULL, zlab = NULL, aspect = c(diff(world.map$range[3:4])/diff(world.map$range[1:2]),
              0.3), panel.aspect = 0.75, lwd = 2, screen = list(z = 30,
              x = -60), par.settings = list(axis.line = list(col = "transparent"),
                box.3d = list(col = "transparent", alpha = 0)))
     

    print(p2)


    # Over US map
    library("maps")
    state.map <- map("state", plot = FALSE, fill = FALSE)

    require(lattice)
    require(latticeExtra)


     # data
     state.info <- data.frame(name = state.name, long = state.center$x,
          lat = state.center$y)


    set.seed(123)
    state.info$yvar<- rnorm (nrow (state.info), 20, 5)


    panel.3dmap <- function(..., rot.mat, distance, xlim,
         ylim, zlim, xlim.scaled, ylim.scaled, zlim.scaled) {
           scaled.val <- function(x, original, scaled) {
          scaled[1] + (x - original[1]) * diff(scaled)/diff(original)
         }
           m <- ltransform3dto3d(rbind(scaled.val(state.map$x,
               xlim, xlim.scaled), scaled.val(state.map$y, ylim,
              ylim.scaled), zlim.scaled[1]), rot.mat, distance)
            panel.lines(m[1, ], m[2, ], col = "grey40")
          }


    pl <- cloud(yvar ~ long + lat, state.info, subset = !(name %in%
           c("Alaska", "Hawaii")), panel.3d.cloud = function(...) {
             panel.3dmap(...)
              panel.3dscatter(...)
     }, col = "blue2",  type = "h", scales = list(draw = FALSE), zoom = 1.1,
                xlim = state.map$range[1:2], ylim = state.map$range[3:4],
              xlab = NULL, ylab = NULL, zlab = NULL, aspect = c(diff(state.map$range[3:4])/diff(state.map$range[1:2]),
              0.3), panel.aspect = 0.75, lwd = 2, screen = list(z = 30,
              x = -60), par.settings = list(axis.line = list(col = "transparent"),
                box.3d = list(col = "transparent", alpha = 0)))
     print(pl)





    Read more »
  • RG#106: add satellite imagery, or open street maps to your plots using openmap package (bing, mapquest)

    library(OpenStreetMap)
    library(rgdal)

    # get world map
    map <- openmap(c(70,-179), c(-70,179))
    plot(map)

    bingmap <- openmap(c(70,-179), c(-70,179), type = "bing")
    plot(bingmap)


    # types available 

    # type = c("osm", "osm-bw", "maptoolkit-topo", "waze", "mapquest", "mapquest-aerial", "bing", "stamen-toner", "stamen-terrain", "stamen-watercolor", "osm-german", "osm-wanderreitkarte", "mapbox", "esri", "esri-topo", "nps", "apple-iphoto", "skobbler", "cloudmade-<id>", "hillshade", "opencyclemap", "osm-transport", "osm-public-transport", "osm-bbike", "osm-bbike-german")


    #zoom maps, plot a portion
    upperLeft, lowerRight
    lat <- c(43.834526782236814, 30.334953881988564)
    lon <- c(-85.8857421875, -70.0888671875)
    southest <- openmap(c(lat[1],lon[1]),c(lat[2],lon[2]),zoom=7,'osm')
    plot(southest) 




    library(UScensus2000tract)
    data(south_carolina.tract)
    lat <- c(35.834526782236814, 30.334953881988564)
    lon <- c(-85.8857421875, -70.0888671875)
    southest <- openmap(c(lat[1],lon[1]),c(lat[2],lon[2]),zoom=7,'osm')
    south_carolina.tract <- spTransform(south_carolina.tract,osm())

    plot(southest)
    plot(south_carolina.tract,add=TRUE,col=(south_carolina.tract@data$med.age>35)+4)


    plot(southest)
    plot(south_carolina.tract,add=TRUE,col=(south_carolina.tract@data$med.age>55)+4)






    Read more »
  • RG#105: Line range or Crossbar plot

    set.seed (1234)
    require(ggplot2)
    Xv <- c(rnorm (500, 10,3), rnorm (500, 50, 20), rnorm (500, 70, 20))
    Yv <- c(rnorm (500, 10,3), rnorm (500, 70, 5), rnorm (500, 30, 5))
    myd <- data.frame (Xv, Yv )

    m <- ggplot(myd, aes(x = Xv, y = Yv)) +
      geom_point() + geom_density2d() + theme_bw()
      
    m


    myd <- data.frame (Xv = c("A", "B", "C", "E", "F"), Yv = c( 35, 35, 69, 15, 60),
                        errs = c(5, 10, 3, 5, 12))

    se <- ggplot(myd, aes(Xv, Yv,
      ymin = Yv - errs, ymax= Yv + errs, colour = Xv))
      
    # line range plot 

    se + geom_linerange() + theme_bw()





    se + geom_linerange()  + coord_flip() + theme_bw()

    se + geom_crossbar(width = 0.5) + theme_bw()



    Read more »
  • RG#104: 2d density plots

    require(ggplot2)

    set.seed (1234)
    Xv <- c(rnorm (500, 10,3), rnorm (500, 50, 20), rnorm (500, 70, 20))
    Yv <- c(rnorm (500, 10,3), rnorm (500, 70, 5), rnorm (500, 30, 5))
    myd <- data.frame (Xv, Yv )

    m <- ggplot(myd, aes(x = Xv, y = Yv)) +
      geom_point() + geom_density2d() + theme_bw()
      
    m



    Read more »
  • RG#103: Combing different types of plot in trellis type

    require(lattice)
    require(latticeExtra)
    #xy plot, conditioned from quakes data, with histogram 
    qk <- xyplot(lat ~ long | equal.count(depth, 2), quakes,
    aspect = "iso", pch = "+", cex = 1.5, xlab = NULL, ylab = NULL)
    qk <- c(qk, depth = histogram(quakes$depth), layout = c(3, 1))
    update(qk, scales = list(at = list(NA, NA, NA, NA),
    y = list(draw = FALSE)))
    # suppress scales on the first 2 panels
    update(qk, scales = list(at = list(NULL, NULL, NA), y = list(draw = FALSE)))
    Read more »
  • RG#102: Double Y axis trellis plot (weather data example)

    require(latticeExtra)
    require (lattice)

    data(SeatacWeather)
    tempatures <- xyplot(min.temp + max.temp ~ day | month,
    data = SeatacWeather, type = "l", layout = c(3, 1))
    rainfall <- xyplot(precip ~ day | month, data = SeatacWeather, type = "h", lwd = 4)

    doubleYScale(tempatures, rainfall, style1 = 0, style2 = 3, add.ylab2 = TRUE,
    text = c("min. T", "max. T", "rain"), columns = 3)


    Read more »
  • RG#101: Plot country map with annotation of regions

    ## Load required packages
    library(maptools)
    library(raster)

     # Download data from gadm.org using getData function 
    admNPL <- getData('GADM', country='NPL', level=3)
    plot(admNPL, bg="lightgreen", axes=T)

    ##Plot 
    plot(admNPL, lwd=10, border="skyblue", add=T)
    plot(admNPL,col="yellow", add=T)
    grid()
    box()

    invisible(text(getSpPPolygonsLabptSlots(admNPL), labels=as.character(admNPL$NAME_3), cex=0.4, col="black", font=2))
    mtext(side=3, line=1, "District Map of Nepal", cex=2)
    mtext(side=1, "Longitude", line=2.5, cex=1.1)
    mtext(side=2, "Latitude", line=2.5, cex=1.1)


    Read more »
  • RG#100: Trellis map plot with heatmap colors

    require(maps)
    library(mapproj)
    worldmap <- map('world', plot = FALSE, fill = FALSE,  projection = "azequalarea")
    country = worldmap$names
    set.seed(1234)
    var.2010 = rnorm (length (country), 20, 10)
    var.2011 = var.2010*1.1 + rnorm (length (country), 0, 1)
    var.2012 = var.2011*0.98 + rnorm (length (country), 0, 4)
    var.2013 = var.2011*0.98 + rnorm (length (country), 0, 30)
    worldt <- data.frame (country, var.2010, var.2011, var.2012, var.2013)
    mapplot(country ~ var.2011, worldt, map = map("world",     plot = FALSE, fill = TRUE))
    mapplot(country ~ var.2010 + var.2011 + var.2012 + var.2013, data = worldt, map = map("world",     plot = FALSE, fill = TRUE))
    # trellis plot for country maps not available in maps package:

     require(maptools)
    # get the map; may need sometime to be loaded 
    con <- url("http://gadm.org/data/rda/NPL_adm3.RData")
    print(load(con))
    close(con)


    # from your data file working directory 
    ## load ("NPL_adm3.RData")

    # data 
    districts = gadm$NAME_3
    set.seed(1234)
    var1 <- rnorm (length (districts), 100, 30)
    var2 <- rnorm (length (districts), 100, 30)
     myd <- data.frame (districts, var1, var2)    







    # US county level map 
    uscountymap <- map('county', plot = FALSE, fill = FALSE,  projection = "azequalarea")
    county = uscountymap$names
    set.seed(1234)
    var.2010 = rnorm (length (county), 50, 10)
    var.2011 = var.2010*1.1 + rnorm (length (county), 0, 5)
    var.2012 = var.2011*0.98 + rnorm (length (county), 0, 10)
    var.2013 = var.2011*1.2 + rnorm (length (county), 0, 15)
    uscounty <- data.frame (county, var.2010, var.2011, var.2012, var.2013)
    mapplot(county ~ var.2010 + var.2011 + var.2012 + var.2013, data = uscounty, map = map("county",    plot = FALSE, fill = TRUE))


    # US state level map 

    usstmap <- map('state', plot = FALSE, fill = FALSE,  projection = "azequalarea")
    state = usstmap$names
    set.seed(1234)
    var.2010 = rnorm (length (state), 50, 10)
    var.2011 = var.2010*1.1 + rnorm (length (state), 0, 5)
    var.2012 = var.2011*0.98 + rnorm (length (state), 0, 10)
    var.2013 = var.2011*1.2 + rnorm (length (state), 0, 15)
    usst <- data.frame (county, var.2010, var.2011, var.2012, var.2013)
    mapplot(state ~ var.2010 + var.2011 + var.2012 + var.2013, data = usst, map = map("state",    plot = FALSE, fill = TRUE), colramp = colorRampPalette(c("green", "purple")))

    Read more »
  • RG#99: cloud 3D bars with heatmap

    require(lattice)
    require(latticeExtra)
    data(VADeaths)
    cloud(VADeaths, panel.3d.cloud = panel.3dbars,
    xbase = 0.4, ybase = 0.4, zlim = c(0, max(VADeaths)),
    scales = list(arrows = FALSE, just = "right"), xlab = NULL, ylab = NULL,
    col.facet = level.colors(VADeaths, at = do.breaks(range(VADeaths), 20),
    col.regions = cm.colors,
    colors = TRUE),
    colorkey = list(col = cm.colors, at = do.breaks(range(VADeaths), 20)),
    screen = list(z = 40, x = -30))
    Read more »
  • RG#98: Horizon plot (time series data)

    require(latticeExtra)
    #example data 
    set.seed(123)
    mydat <- ts(matrix(cumsum(rnorm(150 * 10)), ncol = 10))
    colnames(mydat) <- paste("TS", letters[1:10], sep = "-")

    #simple line plot
    xyplot(mydat, scales = list(y = "same"))


    # panel with different origin and scale:
    horizonplot(mydat, layout = c(1,12), colorkey = TRUE) +  
     layer(panel.scaleArrow(x = 0.99, digits = 1, col = "grey",
                             srt = 90, cex = 0.7)) +
    layer(lim <- current.panel.limits(),
    panel.text(lim$x[1], lim$y[1], round(lim$y[1],1), font = 2,
            cex = 0.7, adj = c(-0.5,-0.5), col = "blue")) 
    Read more »
  • RG#97: Error bar plot with significance (line connecting) - publication purpose

    myd <- data.frame (X = c(1:12,1:12),
                       Y = c(8, 12, 13, 18,  22, 16, 24, 29,  34, 15, 8, 6,
                             9, 10, 12, 18, 26, 28, 28, 30, 20, 10, 9, 9),
                       group = rep (c("A-group", "B-group"), each = 12),
                       error = rep (c(2.5, 3.0), each = 12))

    plot1 = ggplot(data = myd, aes(x=X, y=Y, fill=group, width=0.8) ) +
         geom_errorbar(aes(ymin=Y, ymax=Y+error, width = 0.2), position=position_dodge(width=0.8)) +
         geom_bar(stat="identity", position=position_dodge(width=0.8)) +
         geom_bar(stat="identity", position=position_dodge(width=0.8), colour="black", legend=FALSE) +
         scale_fill_manual(values=c("grey70", "white")) +
         scale_x_discrete("X", limits=c(1:12)) +
         scale_y_continuous("Y (units)", expand=c(0,0), limits = c(0, 40), breaks=seq(0, 40, by=5)) + ggtitle ("My nice plot") +
         theme_bw() +
        theme( plot.title = element_text(face="bold", size=14),
              axis.title.x = element_text(face="bold", size=12),
              axis.title.y = element_text(face="bold", size=12, angle=90),
              panel.grid.major = element_blank(),
              panel.grid.minor = element_blank(),
              axis.text.y=element_text(angle=90, hjust=0.5),
              legend.title = element_blank(),
              legend.position = c(0.85,0.85),
              legend.key.size = unit(1.5, "lines"),
              legend.key = element_rect()
         )

    print(plot1)





       # connecting segments
       plot1 + geom_segment(aes(x=5.7, y=16, xend=5.7, yend=28.5), col = "gray80")+
       geom_segment(aes(x=5.7, y=28.5, xend=6.1, yend=28.5), col = "gray80")      +
       geom_segment(aes(x=6.1, y=28.5, xend=6.1, yend=28), col = "gray80") +
          annotate("text", x=5.7, y=31.2, label="***")





    Read more »
  • RG#96: Basic point and line graph with error bars (publication purpose)

    myd <- data.frame (X = c(1:12,1:12),
                       Y = c(8, 12, 13, 18,  22, 16, 24, 29,  34, 15, 8, 6,
                             9, 10, 12, 18, 26, 28, 28, 30, 20, 10, 9, 9),
                       group = rep (c("A-group", "B-group"), each = 12),
                       error = rep (c(2.5, 3.0), each = 12))
                       
    require(ggplot2)
    require(grid)
    # line and point plot
    f1 = ggplot(data = myd, aes(x = X, y = Y, group = group) )  # lesion becomes a classifying factor
    f2 <- f1 + geom_errorbar(aes(ymin = Y - error, ymax = Y + error), width=0.3) +
    geom_line() + geom_point(aes(shape=group, fill=group), size=5)

     f3 <- f2 +  scale_x_continuous("X (units)", breaks=1:12) +
         scale_y_continuous("Y (units)", limits = c(0, 40), breaks=seq(0, 40, by = 5)) +
         scale_shape_manual(values=c(24,21)) +
         scale_fill_manual(values=c("white","black")) +
         stat_abline(intercept=0, slope=0, linetype="dotted") +
         annotate("text", x=11, y=10, label="X") +
         theme_bw()

       optns <- theme (
              plot.title = element_text(face="bold", size=14),
              axis.title.x = element_text(face="bold", size=12),
              axis.title.y = element_text(face="bold", size=12, angle=90),
              panel.grid.major = element_blank(),
              panel.grid.minor = element_blank(),
              legend.position = c(0.2,0.8),
              legend.title = element_blank(),
              legend.text = element_text(size=12),
              legend.key.size = unit(1.5, "lines"),
              legend.key = element_blank()
         )
    f3 +  ggtitle ( "MY awsome plot for publication") + optns



    Read more »
  • RG#95: Interactive Biplot

    # data
    set.seed(1234)
    P <- vector()
    DF <- as.data.frame(matrix(rep(NA, 100), nrow=10))
    names(DF) <- c(paste("M",1:10, sep=""))
    for(i in 1:10) {
    DF[,i] <- rnorm(10, 10, 3)
    }
    rownames (DF) <- paste("O", 1:10, sep = "")

    require(BiplotGUI)


    Biplots(Data = DF, groups = rep(1, nrow(DF)),PointLabels = rownames(DF),
    AxisLabels = colnames(DF))



    # you can work in the interactive menu window



    Read more »
  • RG#94: Hive plot (diagram)
    require(HiveR)

    # test random data, nodes = 200, edge 200   
    td1 <- ranHiveData(nx = 3, nn = 200, ne = 200)
    plotHive(td1, bkgnd = "cornsilk"
     
     

     
    Read more »
  • RG#94: Arch Diagram
    # 'devtools' need to be installed if not installed 
    install.packages("devtools")

    # load devtools
    library(devtools)

    # install 'arcdiagram'
    install_github('arcdiagram', username='gastonstat')
     
     
     library(arcdiagram)

    # star graph with 20 nodes
    star_graph = graph.star(20, mode="out")

    # extract edgelist
    staredges = get.edgelist(star_graph)

    # default arc diagram
    arcplot(staredges)
      
     
     # different ordering, arc widths, arc colors, and node sizes
    set.seed(1234)
    orderv <- sample (1:20)
    labelv <- paste("nd",1:20,sep="-")
    archlinewd <- 4*runif(20,.5,2)
    colarcs <- heat.colors (19)
    nodesize <- runif(20,1,3)
    arcplot(staredges, ordering=orderv, labels= labelv,
    lwd.arcs= archlinewd, col.arcs= colarcs,
    show.nodes=TRUE, pch.nodes=21, cex.nodes= nodesize,
    col.nodes="lightseagreen", bg.nodes="midnightblue", lwd.nodes= 2)
     
     
     
     
     
     
     
     
    Read more »
  • RG#93: Add countour or heat map plot to XY scatter plot
    # data
    set.seed(1234)
    n <- 10000
    X = rnorm (n, 10, 4)
    Y = X*1.5 + rnorm (n, 0, 8)


    ## colour brewing
    library(RColorBrewer)
    g = 11
    my.cols <- rev(brewer.pal(g, "RdYlBu"))

    #compute 2D kernel density

    # kernel density using MASS 
    library(MASS)
    z <- kde2d(X, Y, n=50)

    plot(X, Y, xlab="X", ylab="Y", pch=19, cex=.3, col = "gray60")
    contour(z, drawlabels=FALSE, nlevels=g, col=my.cols, add=TRUE, lwd = 2)
    abline(h=mean(Y), v=mean(X), lwd=2, col = "black")
    legend("topleft", paste("r=", round(cor(X, Y),2)), bty="n")

    ## estimate the z counts
    prob <- c(.99, .95, .90, .8, .5, .1, 0.05)
    dx <- diff(z$x[1:2])
    dy <- diff(z$y[1:2])
    sz <- sort(z$z)

    c1 <- cumsum(sz) * dx * dy

    levels <- sapply(prob, function(x) {
                  approx(c1, sz, xout = 1 - x)$y })

    plot(X,Y, col = "gray80", pch = 19, cex = 0.3)
    contour(z, levels= round (levels,7), add=T, col = "red")




    # smooth scatter
    require(KernSmooth)
    smoothScatter(X, Y, nrpoints=.3*n, colramp=colorRampPalette(my.cols), pch=19, cex=.3, col = "green1")





     
    Read more »
  • RG#91: Plot bar or pie chart over world map using rworldmap package
    require (rworldmap)
    require(rworldxtra)

    # get world map 

    plot(getMap(resolution = "high" ))


    ##getting example data
    dataf <- getMap()@data 
    mapBars( dataf, nameX="LON", nameY="LAT" , nameZs=c('GDP_MD_EST',
    'GDP_MD_EST','GDP_MD_EST') , mapRegion='asia' , symbolSize=2  ,
     barOrient = 'horiz' )

    mapBars( dataf, nameX="LON", nameY="LAT" , nameZs=c('GDP_MD_EST','GDP_MD_EST',
    'GDP_MD_EST') , mapRegion='asia' , symbolSize=3  , barOrient = 'vert' ,  
    oceanCol = "blue1", landCol = "lightgreen")



     
    mapPies( dataf,nameX="LON", nameY="LAT", nameZs=c('GDP_MD_EST','GDP_MD_EST',
    'GDP_MD_EST','GDP_MD_EST'),mapRegion='asia', oceanCol = "lightseagreen",
     landCol = "gray50")




    Read more »
  • RG#90: fluctutation diagram: graphical representation of a contingency table
    set.seed (1234)
    myd <- data.frame (x1 = rnorm (1000, 15, 5), x2 = sample (c("A", "B", "C"), 1000, replace = TRUE), x3 = sample (c(1,2,2), 1000, replace = TRUE))


    # fluctuation plot
    require (ggplot2)

    ggfluctuation(table(myd$x2, myd$x3))+theme_bw()
    ggfluctuation(table(myd$x2, myd$x3), type="colour")+theme_bw()
    ggfluctuation(table(cut (myd$x1,4), myd$x3), type="colour")+theme_bw()



    ggfluctuation(table(cut (myd$x1,5), myd$x3))+theme_bw()




    Read more »
  • RG#89: Raster in R

     require(grid)
    grid.raster(matrix(colors()[1:600], ncol=20))





    set.seed(1234)

    mt <- matrix (sample(c("red","red1", "yellow", "purple", "green1", "green4", "blue"), 10000, replace = TRUE), ncol = 100)

    grid.raster(mt)


     
    rgb.palette <- colorRampPalette(c("red", "orange", "blue"),
    space = "rgb")
     
    bg.palette <- colorRampPalette(c("blue", "green"), space = "rgb")

    gr.palette <-  colorRampPalette(c("green", "red"),
    space = "rgb")
     
     
    colrs <- matrix(c(rgb.palette(20),bg.palette(20),gr.palette(20)),   nrow=6, ncol=10)
    grid.newpage()
     grid.raster(colrs)
      grid.raster(colrs, interpolate=FALSE)
    # raster in ggplot2
    require(ggplot2)

    # Generate data
    funp <- function (n,r=2) {
    xv <- seq(-r*pi, r*pi, len=n)
    df1 <- expand.grid(x=xv, y=xv)
    df1$r <- sqrt(df1$x^2 + df1$y^2)
    df1$z <- cos(df1$r^2)*exp(-df1$r/6)
    df1
    }
    qplot(x, y, data = funp(1000), fill = z, geom = "raster") + theme_bw()








    Read more »
  • RG#88: Plot pie over google map
    require(ggplot2)
    library(ggmap)
    require(grid)

     # blank theme function
    blk_theme <- function (){theme(
         line = element_blank(),
         rect = element_blank(),
         text = element_blank(),
         axis.ticks.length = unit(0, "cm"),
         axis.ticks.margin = unit(0.01, "cm"),
         legend.position = "none",
         panel.margin = unit(0, "lines"),
         plot.margin = unit(c(0, 0, -.5, -.5), "lines"),
         complete = TRUE
       )
       }
      
    dha = get_map(location = c(lon = 80.56410278, lat = 28.7089375), zoom = 14, maptype = 'roadmap', source = "google")
    dhamp <-  ggmap(dha) + blk_theme ()


    # data 1
    d1 <- sample ( c("A", "B", "C", "D"), 100, replace = TRUE)
    d2 <- sample ( c("A", "B", "C"), 100, replace = TRUE)
    d3 <- sample ( c("A", "B"), 100, replace = TRUE)
    d4 <- sample ( c("A", "B"), 100, replace = TRUE)
    d5 <- sample ( c("A", "B", "C", "D"), 100, replace = TRUE)

    myd <- data.frame(x = factor(1),d1, d2, d3, d4, d5)

    # pie charts
    pie1 <- qplot(x, d1, data = myd, geom = 'bar', fill = d1) +
      coord_polar(theta = 'y') + blk_theme()
    pie2 <- qplot(x, d2, data = myd, geom = 'bar', fill = d2) +
      coord_polar(theta = 'y') + blk_theme()
    pie3 <- qplot(x, d3, data = myd, geom = 'bar', fill = d3) +
      coord_polar(theta = 'y') + blk_theme()
    pie4 <- qplot(x, d4, data = myd, geom = 'bar', fill = d4) +
      coord_polar(theta = 'y') + blk_theme()
    pie5 <- qplot(x, d5, data = myd, geom = 'bar', fill = d5) +
      coord_polar(theta = 'y') + blk_theme()

    # plotting and viewports
    # function
    vplayout <- function(x, y)viewport(layout.pos.row = x, layout.pos.col = y)
    ####
    grid.newpage()
    pushViewport(viewport(layout = grid.layout(20,20)))
    # print mainmap
    print(dhamp, vp = vplayout(1:20, 1:20))

      # need to find manually, maximum and minimum lon and lat
     # click on google earth to find
     maximum.lon <- 80.591
     minimum.lon <- 80.538
     maximum.lat <- 28.735
     minimum.lat <- 28.685

     # X and Y cordinates where pie need to draw
     pieposX <- c(28.72, 28.73, 28.70,28.73,28.69 )
     pieposY <- c(80.55, 80.58, 80.58, 80.55,80.58)

    ycoord <- 20 * ((pieposX - minimum.lat)/(maximum.lat - minimum.lat))
    xcoord <- 20 * ((pieposY - minimum.lon)/(maximum.lon - minimum.lon))

    # plot pie over layer
    print(pie1, vp = vplayout(xcoord[1], ycoord[1]))
    print(pie2, vp = vplayout(xcoord[2], ycoord[2]))
    print(pie3, vp = vplayout(xcoord[3], ycoord[3]))
    print(pie4, vp = vplayout(xcoord[4], ycoord[4]))
    print(pie5, vp = vplayout(xcoord[5], ycoord[5]))



    # bigger pie plot  

    grid.newpage()
    pushViewport(viewport(layout = grid.layout(7,7)))
    # print mainmap
    print(dhamp, vp = vplayout(1:7, 1:7))

      # need to find manually, maximum and minimum lon and lat
     # click on google earth to find
     maximum.lon <- 80.591
     minimum.lon <- 80.538
     maximum.lat <- 28.735
     minimum.lat <- 28.685

     # X and Y cordinates where pie need to draw
     pieposX <- c(28.72, 28.73, 28.70,28.73,28.69 )
     pieposY <- c(80.55, 80.58, 80.58, 80.55,80.58)

    ycoord <- 7 * ((pieposX - minimum.lat)/(maximum.lat - minimum.lat))
    xcoord <- 7 * ((pieposY - minimum.lon)/(maximum.lon - minimum.lon))

    # plot pie over layer 
    print(pie1, vp = vplayout(xcoord[1], ycoord[1]))
    print(pie2, vp = vplayout(xcoord[2], ycoord[2]))
    print(pie3, vp = vplayout(xcoord[3], ycoord[3]))
    print(pie4, vp = vplayout(xcoord[4], ycoord[4]))
    print(pie5, vp = vplayout(xcoord[5], round (ycoord[5],0)))







    Read more »
  • RG#87: histogram / bar chart over map
    library(ggsubplot)
    library(ggplot2)
    library(maps)
    library(plyr)

    #Get world map info
    world_map <- map_data("world")

    #Create a base plot
    p <- ggplot()  + geom_polygon(data=world_map,aes(x=long, y=lat,group=group), col = "blue4", fill = "lightgray") + theme_bw()

    # Calculate the mean longitude and latitude per region (places where subplots are plotted),
    cntr <- ddply(world_map,.(region),summarize,long=mean(long),lat=mean(lat))

    # example data
     myd <- data.frame (region = rep (c("USA","China","USSR","Brazil", "Australia","India", "Canada"),5),
                        categ = rep (c("A", "B", "C", "D", "E"),7), frequency = round (rnorm (35, 8000, 4000), 0))
                       

    subsetcntr  <- subset(cntr, region %in% c("USA","China","USSR","Brazil", "Australia","India", "Canada"))

    simdat <- merge(subsetcntr, myd)
    colnames(simdat) <- c( "region","long","lat", "categ", "myvar" )


     myplot  <- p+geom_subplot2d(aes(long, lat, subplot = geom_bar(aes(x = categ, y = myvar, fill = categ, width=1), position = "identity")), ref = NULL, data = simdat)

    print(myplot)






    Read more »
  • RG#86: 3D XY plot with sphare plots (interactive)
    # data
    myd <-data.frame(name =c("A","B","C","D", "E"),
    var_x=c(6,7,11,1,8),
    var_y=c(9,2,9,4, 2),
    var_z=c(4,1,6,5,1),
    point_size=c(6,3,6,3, 5)
    )
    myd$pradius <- myd$point_size*0.15

    require(rgl)

    spheres3d(myd[,2:4], radius = myd[,6], col = c("darkred", "green", "yellow", "orange", "purple"), alpha = 0.8)
    axes3d(box = TRUE)

    #title3d(xlab = "var_x", ylab = "var_y", zlab = "var_z")
    #text3d(myd[1,2:5], texts = "A")

    segs <- rbind(myd[1:2,2:5], myd[2:3,2:5], myd[3:4,2:5], myd[4:5,2:5], myd[c(5,1),2:5])
    segments3d(segs, col = "blue", lwd = 2)

    # take a SNPshot
    rgl.snapshot ("my3dplot.png", fmt = "png")





    # you can rotate the axis and take another snapshot, and shave as different name
    rgl.snapshot ("my3dplot2.png", fmt = "png")

    Read more »

Copyright Use-R.com 2012 - 2016 ©