#================================================================ # ASSOCIATION NAVIGATOR -- 2013/11/26 # Code by A. Buja, http://stat.wharton.upenn.edu/~buja/ #================================================================ # 2013/06/26 After experimenting on a Mac w Natalia: # .Platform$OS.type # From help(): character string, giving the Operating System (family) of the computer. One of "unix" or "windows". # # 2013/05/21 London, after visiting Catherine in Dublin: # She managed to get the association navigator going under Cairo. # Below are her comments. Might still be worth contacting # Michael Lawrence and ask about the slowness. ## yes, the cairoDevice works on the mac. Tried it with the data ## below which has 1600 rows (used as variables), ## not nearly as smooth I think as the windows version though. ## The problem last night was that your code was executing dev.new() ## which made a new device of the default type. ## Once I reset the default to be cairo it worked ok. ## require(cairoDevice) ## source("http://stat.wharton.upenn.edu/~buja/association-navigator.R") ## require(HSAUR2) ## x <- CHFLS ## x <- x[,sapply(x, is.numeric)] ## x<- t(as.matrix(x)) ## options(device=function() Cairo(width=8,height=8)) ## a.n <- a.nav.create(x) ## a.nav.run(a.n) #================================================================ # DEVELOPMENT NOTES: # Fix: z scatterplot matrix bombs when highlight outside of range; prevent highlight outside range # See end of file for more development notes and execution code # To do: Adjustment clean-up: access to 'cors.compl' and 'notna' needs to cohere to adjusted vars. # for uninformative variables (one or zero value only other than NAs), skip adjustment. # Project: Trimmed correlations, trimming, 1, 2, 3, ... cases from each extreme, randomly for ties # Machinery: gaga <- c(NA,NA,1,2,1,2,1,2,NA,NA); gaga[rank(gaga, ties="random", na.last="keep")==6 & !is.na(gaga)] # Project: replace correlations of data with correlations of correlations against a selected set of variables (e.g., SRS) # Project: Bayes net #=============================== ASSOCIATION NAVIGATOR === R-CODE ================================ #=============================== Begin Experimental =============================== # Selector/creator of ANs: ## # Find numeric matrices and create a factor that tells whether they have same nrow() and rownames(): ## tmp <- sapply(ls(), function(nam) { is.matrix(get(nam)) & is.numeric(get(nam)) }) ## tmp.nams <- names(tmp)[tmp] ## tmp.rn <- sapply(tmp.nams, function(nam) { paste(paste(c(nrow(get(nam)),rownames(get(nam))), collapse="|"),"|",sep="") }) ## substr(tmp.rn, 1, 20) # Excellent! Matrices w/o rownames have only nrow ## as.numeric(factor(tmp.rn)) # Factor that groups matrices matched by nrow() and rownames() ## # Want them in descending size of group: ## tmp.fac <- as.numeric(factor(tmp.rn)) ## # Reverse sort by size of levels: ## tmp.rk <- max(tmp.fac)+1-rank(table(tmp.fac), ties="first") ## tmp.fac <- as.vector(tmp.rk[tmp.fac]) ## table(tmp.fac) ## ord <- order(tmp.fac) ## tmp.fac[ord] ## data.frame(nrow.and.init.rownames=substr(tmp.rn[ord],1,20), grp=tmp.fac[ord]) ## sel.lst <- paste(format(names(tmp.rn)), ## format(sub("[|].*","",tmp.rn)), ## format(sapply(names(tmp.rn), function(nam) ncol(get(nam)))), ## format(substr(sub("[^|]*[|]","",tmp.rn[ord]),1,30)), ## format(tmp.fac) ## )[ord] ## select.list(sel.lst, multiple=T, title="List of Numeric Matrices") #=============================== End Experimental =============================== # Create an 'association navigator': a.nav.create <- function(dat, remove=FALSE, bl.sep="_") { # Separator for postfix: store locally bl.sep <- bl.sep # For restoring later. options.orig <- options() # Coerce to matrix in case 'dat' is a dataframe. dat <- as.matrix(dat) if(!is.numeric(dat)) { cat("a.nav.create(): Data matrix is not numeric."); return() } # Fix colnames if necessary: if(is.null(colnames(dat))) { # There are no colnames; create them: colnames(dat) <- paste("Var", 1:ncol(dat)) } else if(length(table(table(colnames(dat)))) != 1) { # Make colnames unique: colnames(dat) <- paste("(", 1:ncol(dat), ") ", colnames(dat), sep="") } # Fix rownames if necessary: if(is.null(rownames(dat))) { # There are no rownames; create them: rownames(dat) <- paste("Obs", 1:nrow(dat)) } else if(length(table(table(rownames(dat)))) != 1) { # Make colnames unique: rownames(dat) <- paste("(", 1:nrow(dat), ") ", rownames(dat), sep="") } # Remove uninformative variables if(remove) { cat("a.nav.create(): Removing uninformative variables (all NA or constant)... \n") sds <- apply(dat, 2, sd, na.rm=T) dat <- dat[,!is.na(sds) & sds!=0, drop=F] if(ncol(dat) <= 1) { cat("a.nav.create(): No informative variables\n"); return() } } Nrow <- nrow(dat) Ncol <- ncol(dat) dat.orig <- dat # backup data because we might do trimming on 'dat' cat("Calculating joint NA patterns... \n") notna <- crossprod(!is.na(dat)) # number of complete pairs # Needs regeneration after each random ranking npairs <- 3 # minimum number of pairs for correlations to be shown ntrim <- 0 # number of cases trimmed from both extremes of all variables ## cat("Calculating mean-imputed correlations... \n") ## std <- scale(dat); std[is.na(std)] <- 0 ## sds <- sqrt(colSums(!is.na(dat))); sds[sds==0] <- NA ## cors.imp <- cor(std)*sds*rep(sds,rep(Ncol,Ncol))/notna; rm(std) # mean-imputed cor; ## cors.imp[notna<=3] <- NA ## # The above 'cors.imp' are 'shrunk'; properly standardized they would be cors.imp*(Nrow-1)/(notna-1) ## # but this normalization does not work; it blows up for small 'notna'. cat("Calculating correlations from complete pairs... \n") options(warn=-1) cors.compl <- cor(dat, use="pairwise.complete.obs") # Expensive and NA-prone under NULL options(warn=1) cors.compl[notna 1) { blot.hist.cur - 1 } else { blot.hist.sz } blot.xlim <<- blot.hist.xlim[,blot.hist.cur] blot.ylim <<- blot.hist.ylim[,blot.hist.cur] }, # back in blockplot window lim history "W" = { blot.hist.cur <<- if(blot.hist.cur < blot.hist.sz) { blot.hist.cur + 1 } else { 1 } blot.xlim <<- blot.hist.xlim[,blot.hist.cur] blot.ylim <<- blot.hist.ylim[,blot.hist.cur] }, # back in blockplot window lim history "H" = { tmp <<- winDialogString("Enter search string/regexpr to find a horizontal variable:", default=strips.regexp) # Can return NULL if(!is.null(tmp)) { strips.regexp <<- tmp var.names <- grep(strips.regexp, colnames(dat), val=T) if(length(var.names)>0) { if(length(var.names)>1) { # select.list() returns "" on cancel or ESC, else the selected string var.name <- select.list(var.names, multiple=F, title="Choose a Horizontal Variable:") } else { var.name <- var.names } # length(var.names)==1 if(nchar(var.name)>0) { # Need to check again because select.list() above could return "" i.var <- match(var.name, colnames(dat)) # no match: returns NA, else the first match index if(!is.na(i.var)) { # Should not be necessary, but still match() could return NA blot.xlim <<- blot.xlim - mean(blot.xlim) + i.var a.nav.blot.hist.incr() strips.x <<- c(strips.x, i.var) xhair.xy[1] <<- i.var } } } } }, "V" = { tmp <<- winDialogString("Enter search string/regexpr to find a vertical variable:", default=strips.regexp) # Can return NULL if(!is.null(tmp)) { strips.regexp <<- tmp var.names <- grep(strips.regexp, colnames(dat), val=T) if(length(var.names)>0) { if(length(var.names)>1) { # select.list() returns "" on cancel or ESC, else the selected string var.name <- select.list(var.names, multiple=F, title="Choose a Vertical Variable:") } else { var.name <- var.names } # length(var.names)==1 if(nchar(var.name)>0) { # Need to check again because select.list() above could return "" i.var <- match(var.name, colnames(dat)) # no match: returns NA, else the first match index if(!is.na(i.var)) { # Should not be necessary, but still match() could return NA blot.ylim <<- blot.ylim - mean(blot.ylim) + i.var strips.y <<- c(strips.y, i.var) xhair.xy[2] <<- i.var } } } } }, "@" = { tmp <<- winDialogString("Enter search string/regexpr to find a variable, hor&vert:", default=strips.regexp) # Can return NULL if(!is.null(tmp)) { strips.regexp <<- tmp var.names <- grep(strips.regexp, colnames(dat), val=T) if(length(var.names)>0) { if(length(var.names)>1) { # select.list() returns "" on cancel or ESC, else the selected string var.name <- select.list(var.names, multiple=F, title="Choose a Variable:") } else { var.name <- var.names } # length(var.names)==1 if(nchar(var.name)>0) { # Need to check again because select.list() above could return "" i.var <- match(var.name, colnames(dat)) # no match: returns NA, else the first match index if(!is.na(i.var)) { # Should not be necessary, but still match() could return NA blot.xlim <<- blot.xlim - mean(blot.xlim) + i.var blot.ylim <<- blot.ylim - mean(blot.ylim) + i.var a.nav.blot.hist.incr() strips.x <<- c(strips.x, i.var) strips.y <<- c(strips.y, i.var) xhair.xy[1:2] <<- i.var } } } } }, "ctrl-C" = { rect.xxyy <<- matrix(NA, ncol=4, nrow=0, # clean rectangles dimnames=list(NULL,c("x1","x2","y1","y2"))) }, "ctrl-D" = { strips.x <<- c(); strips.y <<- c() }, # clean rectangles "ctrl-X" = { blot.mar <<- c(8,8,.5,.5) # reset blockplot blot.xlim <<- c(0,nrow(tab)+1) blot.ylim <<- c(0,ncol(tab)+1) a.nav.blot.hist.incr() }, "ctrl-Z" = { blot.mar <<- c(8,8,.5,.5); # reset blockplot zoom/pan preserving square aspect ratio blot.xlim <<- c(0,nrow(tab)+1) * max(1,par()$pin[1]/par()$pin[2]) blot.ylim <<- c(0,ncol(tab)+1) * max(1,par()$pin[2]/par()$pin[1]) a.nav.blot.hist.incr() }, "0" = { ref.nx <- if(ref.flag) length(ref.vars.x) else 0 ref.ny <- if(ref.flag) length(ref.vars.y) else 0 if(xhair.xy[1]>=blot.xlim[1]+ref.nx & xhair.xy[2]>=blot.ylim[1]+ref.ny) { is <- a.nav.in.rects(xhair.xy, rect.xxyy) if(length(is)>0) rect.xxyy <<- rect.xxyy[-is,,drop=F] } }, "1" = { ref.nx <- if(ref.flag) length(ref.vars.x) else 0 ref.ny <- if(ref.flag) length(ref.vars.y) else 0 if(xhair.xy[1]>=blot.xlim[1]+ref.nx & xhair.xy[2]>=blot.ylim[1]+ref.ny) { rect.xxyy <<- rbind(rect.xxyy, round(xhair.xy)[c(1,1,2,2)]) # Start new focus rectangle i <- nrow(rect.xxyy) # We know i>0 rect.xxyy[i,] <<- pmax(pmin(rect.xxyy[i,], ncol(dat)), 1) } }, "2" = { i <- nrow(rect.xxyy) # Place opposite corner of new focus rect ref.nx <- if(ref.flag) length(ref.vars.x) else 0 ref.ny <- if(ref.flag) length(ref.vars.y) else 0 if(i>0 & xhair.xy[1]>=blot.xlim[1]+ref.nx & xhair.xy[2]>=blot.ylim[1]+ref.ny) { rect.xxyy[i,c(2,4)] <<- round(xhair.xy) # (Not assuming sorted xxyy) rect.xxyy[i,] <<- pmax(pmin(rect.xxyy[i,], ncol(dat)), 1) } }, "3" = { a.nav.focus.rect() }, # focus rectangle: variable block combination under x-hair "4" = { i <- nrow(rect.xxyy) # snap to last focus rectangle maintaining square aspect ratio if(i>0) { # doesn't work properly in the presence of reference variables (see also "5") ii <- a.nav.in.rects(xhair.xy, rect.xxyy) # Is(/are) existing rectangle(s) pointed at by mouse? if(length(ii)>0) i <- ii[length(ii)] # Yes; pick last rectangle; else leave i at last rectangle blot.xlim <<- sort(rect.xxyy[i,1:2])+c(-1,+1)*.75 #+c(-length(ref.vars.x,0) blot.ylim <<- sort(rect.xxyy[i,3:4])+c(-1,+1)*.75 #+c(-length(ref.vars.y,0) if(diff(blot.xlim)0) { # doesn't work properly in the presence of reference variables (see also "5") ii <- a.nav.in.rects(xhair.xy, rect.xxyy) # Is(/are) existing rectangle(s) pointed at by mouse? if(length(ii)>0) i <- ii[length(ii)] # Yes; pick last rectangle; else leave i at last rectangle blot.xlim <<- sort(rect.xxyy[i,1:2])+c(-1,+1)*.75; blot.ylim <<- sort(rect.xxyy[i,3:4])+c(-1,+1)*.75 a.nav.blot.hist.incr() } }, "R" = { ref.vars.x <<- strips.x; # Turn marked variables into reference variables ref.vars.y <<- strips.y ref.flag <<- T }, "r" = { if(mode=="Blockplot") { ref.flag <<- !ref.flag # Toggle showing of reference variables } else if(mode=="Scatterplot") { brush <<- c(mouse.xy[1], mouse.xy[1], mouse.xy[2], mouse.xy[2]) # Reset brush size a.nav.scatterplot() } }, "c" = { a.nav.xhair.draw(col="white") # cross-hair size - xhair.sz <<- xhair.sz/2 a.nav.xhair.draw(col="black") }, "C" = { xhair.sz <<- xhair.sz*2 # cross-hair size + a.nav.xhair.draw() }, "PgUp" = { if(lens.i>=1) { lens.i <<- max(lens.i-1, 1) # step to beginning of lenses a.nav.scatterplot() a.nav.xhair.move(c(which(colnames(tab) %in% lens.lst.xy[[lens.i]][1]), which(rownames(tab) %in% lens.lst.xy[[lens.i]][2]) ) ) } }, "PgDn" = { if(lens.i<=length(lens.lst.xy)) { lens.i <<- min(lens.i+1, length(lens.lst.xy)) # step to end of lenses a.nav.scatterplot() a.nav.xhair.move(c(which(colnames(tab) %in% lens.lst.xy[[lens.i]][1]), which(rownames(tab) %in% lens.lst.xy[[lens.i]][2]) ) ) } }, "End" = { if(length(lens.lst.xy)>0) { lens.i <<- length(lens.lst.xy) # to end of lenses a.nav.scatterplot() a.nav.xhair.move(c(which(colnames(tab) %in% lens.lst.xy[[lens.i]][1]), which(rownames(tab) %in% lens.lst.xy[[lens.i]][2]) ) ) } }, "Home" = { if(length(lens.lst.xy)>0) { lens.i <<- 1; # to beginning of lenses a.nav.scatterplot() a.nav.xhair.move(c(which(colnames(tab) %in% lens.lst.xy[[lens.i]][1]), which(rownames(tab) %in% lens.lst.xy[[lens.i]][2]) ) ) } }, "Del" = { if(length(lens.lst.xy)>1 & length(lens.lst.xy)>=lens.i) { # remove a lens plot from history lens.lst.xy <<- lens.lst.xy[-lens.i] lens.i <<- min(lens.i, length(lens.lst.xy)) a.nav.scatterplot() a.nav.xhair.move(c(which(colnames(tab) %in% lens.lst.xy[[lens.i]][1]), which(rownames(tab) %in% lens.lst.xy[[lens.i]][2]) ) ) } }, "s" = { if(mode=="Blockplot") { mode <<- "Scatterplot" # Toggle blockplot <-> scatterplot brush.reset <<- "Reset" a.nav.scatterplot() } else if(mode=="Scatterplot") { mode <<- "Blockplot" a.nav.scatterplot() a.nav.blockplot() } }, "S" = { if(mode=="Scatterplot") # Switch brushing colors; switch to scatterplot if necessary. brush.col <<- brush.cols[1 + (which(brush.cols %in% brush.col) %% length(brush.cols))] mode <<- "Scatterplot" a.nav.scatterplot() }, "x" = { a.nav.lens.add() }, # Show that block in Lens and add to record of examined plots. "l" = { curve.draw <<- !curve.draw; a.nav.scatterplot() }, # Toggle smooth curve/mean trace "y" = { lens.hor.vert <<- !lens.hor.vert; a.nav.scatterplot() }, # Toggle x<->y in scatterplots "z" = { a.nav.splom() }, # Scatterplot matrix "Z" = { a.nav.splom(equalaxes=T) }, # Scatterplot matrix equal axes "ctrl-O" = { a.nav.cors.trans() tab.has <<- "Correlations\n(Compl.Pairs)" }, "ctrl-I" = { a.nav.pvals.trans() tab.has <<- "P-values\n(Normal)" }, ">" = { if(length(grep("P-values", tab.has))>0) # Select stronger correlations or smaller p-values { pvals.thresh.i <<- min(pvals.thresh.i+1, length(pvals.thresh.lst)) a.nav.pvals.trans() } else { cors.thresh.i <<- min(cors.thresh.i+1, length(cors.thresh.lst)) a.nav.cors.trans() } }, "<" = { if(length(grep("P-values", tab.has))>0) # Allow weaker correlations or larger p-values { pvals.thresh.i <<- max(pvals.thresh.i-1, 1) a.nav.pvals.trans() } else { cors.thresh.i <<- max(cors.thresh.i-1, 1) a.nav.cors.trans() } }, "ctrl-T" = { str <- "Enter number of extreme cases to be trimmed\nfrom both sides of all variables:" trm <- as.integer(winDialogString(str, default=as.character(ntrim))) if(length(trm)>0 && trm>=0 && trm<=(nrow(dat)/2-1) && !(trm==0 && ntrim==0)) { cat("\nTrimming the data for",trm,"... \n") ntrim <<- trm # 'random' ranks for ties, hence different every time (also preserve NAs with 'keep'): # !!!!!!!!!!!! Not good -- should have reproducible results; problem for categoricals # !!!!!!!!!!!! Write code to leave categories untouched unless unless #in.categ<=ntrim ## dat.rank <- apply(dat.orig, 2, rank, ties="random", na.last="keep") ## dat.rank <- apply(dat.rank, 2, function(x) pmin(x, max(x)+1-x)) ## dat <<- dat.orig ## dat[dat.rank<=ntrim] <<- NA # Trimming! cat("\nCalculating trimming quantiles... ") eps <- 1E-50 dat.quantiles <- apply(dat.orig, 2, quantile, c(ntrim+eps,nrow(dat)+1-ntrim-eps)/nrow(dat), na.rm=T) cat("Trimming by setting to NA... ") dat <<- dat.orig for(i in 1:nrow(dat)) { dat[i, dat[i,]dat.quantiles[2,]] <- NA } dat <<- dat cat("Calculating joint NA patterns... \n") notna <<- crossprod(!is.na(dat)) # Number of complete pairs ## cat("Calculating mean-imputed correlations... \n") ## std <- scale(dat); std[is.na(std)] <- 0 ## sds <- sqrt(colSums(!is.na(dat))); sds[sds==0] <- NA ## cors.imp <<- cor(std)*sds*rep(sds,rep(Ncol,Ncol))/notna ## cors.imp[notna0 && prs>=0 && prs<=Nrow && prs!=npairs) { npairs <<- prs switch(tab.has, "Correlations\n(Compl.Pairs)" = { a.nav.kbd("ctrl-O") }, "P-values\n(Normal)" = { a.nav.kbd("ctrl-I") }, "#Complete\nPairs/N" = { a.nav.kbd("ctrl-N") }, "#Missing\nPairs/N" = { a.nav.kbd("ctrl-J") }, "Correlations\n(All Random)" = { a.nav.kbd("ctrl-R") }, "Correlations\n(Block Random)" = { a.nav.kbd("ctrl-B") } ) } }, "ctrl-A" = { tab.abs <<- !tab.abs }, # Ignore sign (all blue) "ctrl-R" = { dat.sim <- dat # Random corrs, all variables permuted dat.sim <- apply(dat, 2, function(x) { sel <- !is.na(x); if(sum(sel)>1) x[sel] <- sample(x[sel]); x }) options(warn=-1) tab[,] <<- cor(dat.sim, use="pairwise.complete.obs") options(warn=1) tab[,] <<- round(abs(tab)^blot.pow * sign(tab),10); tab[notna1) x[sel] <- sample(x[sel]) x }) dat.sim <- sapply(1:length(postfx), function(i) dat.sim[blocks.perm[,postfx[i]],i]) colnames(dat.sim) <- colnames(dat) options(warn=-1) tab[,] <<- cor(dat.sim, use="pairwise.complete.obs") options(warn=1) tab[,] <<- round(abs(tab)^blot.pow * sign(tab), 10) tab[notna0) " Or select ADJUSTORS from the following list:", adjust.ors.varnames, "__________________________________________________________________________________", "", " 'ADJUSTEES' / DEPENDENT y-VARIABLES TO BE ADJUSTED: Search or select", pick.ees, if(length(adjust.ees.varnames)>0) " Or select ADJUSTEES from the following list:", adjust.ees.varnames, ""), preselect=NULL, # c(pick.ors, pick.ees), multiple=T, # Permits length(picks)==0 (nothing selected) on 'cancel' or ESC title="Adjustment: Independent x-Variables To Adjust For, Dependent y-Variables To Be Adjusted") if(any(nchar(picks)>0)) { print(picks) # Adjustors: leave untouched if neither string/regexp nor selection from list if(any(picks==pick.ors)) { tmp <- winDialogString("Enter a string or regexp for ADJUSTOR variables:", default=adjust.ors.regexp) if(!is.null(tmp)) { adjust.ors.regexp <<- tmp adjust.ors.varnames <<- grep(adjust.ors.regexp, colnames(dat), val=T) adjust.ees.varnames <<- setdiff(adjust.ees.varnames, adjust.ors.varnames) # keep them disjoint } } else if(length(intersect(picks, adjust.ors.varnames))>0) { # Adjustor varnames have been subselected adjust.ors.varnames <<- intersect(picks, adjust.ors.varnames) adjust.ees.varnames <<- setdiff(adjust.ees.varnames, adjust.ors.varnames) # keep them disjoint } # Adjustees: leave untouched if neither string/regexp nor selection from list if(any(picks==pick.ees)) { tmp <- winDialogString("Enter a string or regexp for ADJUSTEE variables:", default=adjust.ees.regexp) if(!is.null(tmp)) { adjust.ees.regexp <<- tmp adjust.ees.varnames <<- grep(adjust.ees.regexp, colnames(dat), val=T) } } else if(length(intersect(picks, adjust.ees.varnames))>0) { # Adjustee varnames have been subselected adjust.ees.varnames <<- intersect(picks, adjust.ees.varnames) } # Do adjustment: if(any(picks==pick.do) & length(adjust.ors.varnames)>0 & length(adjust.ees.varnames)>0) { adjust.dat <<- matrix(NA, nrow=nrow(dat), ncol=length(adjust.ees.varnames), dimnames=list(rownames(dat),adjust.ees.varnames)) X <- dat[,adjust.ors.varnames] cat("Adjusting...\n") for(nam in colnames(adjust.dat)) { y <- dat[,nam] if(sum(!is.na(y+X))>10) { res <- resid(lm(y ~ X)) if(length(res)>0) adjust.dat[names(res),nam] <<- res } cat(nam,": length=",length(res)," #NA=",sum(is.na(res))," has names: ", !is.null(names(res)),"\n",sep="") } cat("Calculating adjusted correlations ...") options(warn=-1) adjust.cors <<- cor(dat, adjust.dat, use="pairwise.complete.obs") options(warn=1) adjust.notna <<- crossprod(!is.na(dat), !is.na(adjust.dat)) adjust.flag <<- T if(length(grep("P-values", tab.has))>0) { a.nav.pvals.trans() } else { a.nav.cors.trans() } cat("Done adjusting correlations ...") break } # Undo adjustment: if(any(picks==pick.undo)) { adjust.flag <<- F if(length(grep("P-values", tab.has))>0) { a.nav.pvals.trans() } else { a.nav.cors.trans() } break } # Mark adjustment variables: if(any(picks==pick.mark)) { strips.x <<- na.omit(match(adjust.ors.varnames, colnames(dat))) # na.omit() should not be necessary strips.y <<- na.omit(match(adjust.ees.varnames, colnames(dat))) break } # Sort Adjustees: no break if(any(picks==pick.sort & length(adjust.ors.varnames)>0 & length(adjust.ees.varnames)>0)) { adjust.ees.RSquared <- c() for(v in adjust.ees.varnames) { adjust.ees.RSquared[v] <- summary(lm(y ~ ., data=data.frame(y=dat[,v], dat[,adjust.ors.varnames])))$r.squared cat(round(adjust.ees.RSquared[v],3)," -- ",v,"\n") } ord <- order(adjust.ees.RSquared, decreasing=T) adjust.ees.RSquared <- adjust.ees.RSquared[ord] adjust.ees.varnames <- adjust.ees.varnames[ord] select.list(paste(format(round(adjust.ees.RSquared,3))," ",format(adjust.ees.varnames)), preselect=NULL,title="Sorted List of dependent Adjustees with RSquared") } ## if(any(picks==pick.quit)) break } else break # if(any(nchar(picks)>0)) { } # repeat{ }, "h" = { m.sel <- select.list(a.nav.help, title="BLOCKPLOT HELP/SELECT: mouse/keyboard operation") if(length(grep(" Hit",m.sel))>0) { # Selection corresponds to a key hit m.strip <- sub(":.*$","",sub("^.*Hit ","",m.sel)) # Wipe out from colon to end and " Hit " prefix if(length(grep(",",m.strip))>0) { # Multiple keys such as 'i,I,o,O'? m.key <- select.list(unlist(strsplit(m.strip,",")), title=m.sel) # Yes; ask which. } else { m.key <- m.strip } # No, it was just one key. if(m.strip=="arrows") { m.key <- select.list(c("Up","Down","Left","Right"), title=m.sel) } if(m.strip=="space") { m.key <- " " } if(m.strip=="Q") { return("Quit") } # Quitting if(m.strip!="h") { a.nav.kbd(m.key) } # Avoid recursion of help. } }, # help "?" = { menu(a.nav.info, graphics=TRUE, title="INFO:") }, # User-supplied info about data "Q" = { options(options.orig); return("Quit") } # Quitting ) if(!(key %in% c("PgUp","PgDn","End","Home","Del","S","s","x","y","l","z","Z","l","P","J","G","?","h"))) { mode <<- "Blockplot"; a.nav.blockplot() } if(exists("dev.flush")) dev.flush() NULL # Not quitting 'getGraphicsEvent()' } # ------------------------ End 'a.nav.kbd()' # a.nav.blot.hist.incr <- function() { # Increase blockplot lim-history blot.hist.cur <<- if(blot.hist.cur < blot.hist.sz) { blot.hist.cur + 1 } else { 1 } blot.hist.xlim[,blot.hist.cur] <<- blot.xlim blot.hist.ylim[,blot.hist.cur] <<- blot.ylim ## cat("a.nav.blot.hist.incr: ", blot.hist.cur, blot.hist.xlim[,blot.hist.cur], blot.hist.ylim[,blot.hist.cur], "\n") } # ------------------------ End 'a.nav.blot.hist.incr()' # a.nav.cors.trans <- function() { # Fills 'tab' with correlations, adjusted if necessary, power transformed if(verbose) cat("a.nav.cors.trans\n") tab[,] <<- cors.compl if(adjust.flag) { tab[,adjust.ees.varnames] <<- adjust.cors tab[adjust.ees.varnames,] <<- t(adjust.cors) } tab[,] <<- abs(tab)^blot.pow * sign(tab) # Corrs from complete pairs, power transform, e.g., blot.pw==.7 thresh <- as.numeric(cors.thresh.lst[cors.thresh.i]) if(thresh>0) tab[abs(tab)0; notna incorrect when 'adjust.flag' is set; need adjust.notna ## if(pval.type=="perm") { ## tab <<- ((pvals+pvals.nsim)/2+1)/(pvals.nsim+1) } else ## if(pval.type=="normal") { ## tab <<- 2*pnorm(abs(cors.imp),m=0,s=sqrt(psqrs/notna),low=F) } else ## tab[notna==0] <- 1 } thresh <- as.numeric(pvals.thresh.lst[pvals.thresh.i]) # initially: 0.05 tab[,] <<- pmax(0, 1 - tab/thresh)^blot.pow * sign(cors.compl) tab[notna0) { if(verbose) cat(button) mouse.button <<- button mouse.xy <<- c(grconvertX(x,"ndc","user"), grconvertY(y,"ndc","user")) # dito mouse.screen.xy <<- c(x,y) if(mode=="Blockplot") { xhair.xy <<- mouse.xy a.nav.blockplot() } else if(mode=="Scatterplot") { # Used to reset brush -- no longer needed? xxxxxxxxxxxx mouse.xy <<- c(grconvertX(x,"ndc","user"), grconvertY(y,"ndc","user")) # dito mouse.button <<- button a.nav.scatterplot() } } else { mouse.button <<- -1 } ## cat("|",button) if(exists("dev.flush")) dev.flush() NULL } # ------------------------ End 'a.nav.mousedown()' # # Explorations for slider bar on sides 3 and 4: a.nav.mousemove <- function(button, x, y) { # panning by dragging if(exists("dev.hold")) dev.hold() if(verbose) cat(".") if(length(button)>0) { if(verbose) cat(button,"|") mouse.button <<- button mouse.xy <<- c(grconvertX(x,"ndc","user"), grconvertY(y,"ndc","user")) # Caution: Must not use difference of mouse here! They stem from previous calls to par()!!! dx <- grconvertX(x,"ndc","user") - grconvertX(mouse.screen.xy[1],"ndc","user") dy <- grconvertY(y,"ndc","user") - grconvertY(mouse.screen.xy[2],"ndc","user") mouse.screen.xy <<- c(x,y) # Belongs here AFTER computation of 'dx', 'dy'!!!! if(mode=="Blockplot") { # Dragging if(button==0) { # Dragging blot.xlim <<- blot.xlim - dx blot.ylim <<- blot.ylim - dy a.nav.blot.hist.incr() xhair.xy <<- mouse.xy } a.nav.blockplot() } else if(mode=="Scatterplot") { a.nav.scatterplot() # Does brushing if on device 1. } } else { mouse.button <<- -1 } ## cat("-",button) # Strange: mousemove keeps firing even when the mouse doesn't move if(exists("dev.flush")) dev.flush() NULL } # ------------------------ End 'a.nav.mousemove()' # a.nav.blockplot <- function(hdcopy=F) { ## cat("a.nav.blockplot: blot.pch =",blot.pch,"\n") if(verbose) cat("\na.nav.blockplot") if(!hdcopy) a.nav.to.main.window() if(ref.flag) { ref.ny <- length(ref.vars.y) ref.nx <- length(ref.vars.x) } else { ref.ny <- 0 ref.nx <- 0 } nx <- 1:nrow(tab) ny <- 1:ncol(tab) nxsel <- which(blot.xlim[1]+.5+ref.nx <= nx & nx <= blot.xlim[2]-.5) nysel <- which(blot.ylim[1]+.5+ref.ny <= ny & ny <= blot.ylim[2]-.5) # Convention: tab[x,y] -- First index 'x' runs fastest. tb <- tab[nxsel, nysel] ## # Original selection version: doesn't work when there are lots of ties at the upper end ## sel <- order(abs(tb), decreasing=T)[1:min(length(tb),max.blocks)] ## # Second version based on subsampling: Excellent performance -- avoid 'order()' !!! tb.len <- length(tb) sel <- which(abs(tb) >= quantile(abs(tb), prob=1-min(1,max.blocks/tb.len), na.rm=T, names=F)) tb <- tb[sel] # Convention: Run x (hor) fast, y (ver) slow xgrid <- rep(nxsel,length(nysel))[sel] ygrid <- rep(nysel,rep(length(nxsel),length(nysel)))[sel] ## xycex <- (abs(tb))^.7/(max(c(length(nxsel),length(nysel)))/blot.sz) # blot.sz <- 90 initially xycex <- abs(tb)/(max(diff(blot.xlim),diff(blot.ylim))/blot.sz) # to check whether power is costly xycol <- if(tab.abs) { blot.col.pos } else { c(blot.col.neg,0,blot.col.pos)[2+sign(tb)] } # Region for non-reference variable blocks: usr <- c(blot.xlim[1]+ref.nx, blot.xlim[2], blot.ylim[1]+ref.ny, blot.ylim[2]) # Does this belong here? Shouldn't it be in another function? ## if(mouse.button==1) { ## 2010/11/14 got rid of middle mouse button ## if(usr[1]<=mouse.xy[1] & mouse.xy[1]<=usr[2]) strips.x <<- a.nav.su(c(strips.x, round(mouse.xy[1]))) ## if(usr[3]<=mouse.xy[2] & mouse.xy[2]<=usr[4]) strips.y <<- a.nav.su(c(strips.y, round(mouse.xy[2]))) ## } else if(mouse.button==2) { ## 2010/11/14 no longer accumulation, now toggling of strips -- 2011/09/11 bad, back to accum. if(usr[1]<=mouse.xy[1] & mouse.xy[1]<=usr[2]) { this.strip.x <- round(mouse.xy[1]) if(this.strip.x %in% strips.x) { # strips.x <<- setdiff(strips.x, this.strip.x) # Remove from strips xxxxxxxx } else { strips.x <<- sort(c(strips.x, this.strip.x)) } } if(usr[3]<=mouse.xy[2] & mouse.xy[2]<=usr[4]) { this.strip.y <- round(mouse.xy[2]) if(this.strip.y %in% strips.y) { # strips.y <<- setdiff(strips.y, this.strip.y) # Remove from strips xxxxxxxxx } else { strips.y <<- sort(c(strips.y, this.strip.y)) } } } ## cat("... before plot() in a.nav.blockplot() ") ## print(dev.list()) plot(x=blot.xlim, y=blot.ylim, type="n", xaxt="n", yaxt="n", xlab="", ylab="", xaxs="i", yaxs="i", xlim=blot.xlim, ylim=blot.ylim, bty="n") ## cat("... after plot() in a.nav.blockplot() ") # Draw diagonal squares to mark the instrument blocks: rect(xleft =pmax(var.grps[1,], usr[1]), ybottom =pmax(var.grps[1,], usr[3]), xright =pmax(pmin(var.grps[2,], usr[2]), usr[1]), ytop =pmax(pmin(var.grps[2,], usr[4]), usr[3]), col =rect.diag.col, border=NA) # Draw the current focus rectangles for a pair of instruments: if(nrow(rect.xxyy)>0) rect(pmax(pmin(rect.xxyy[,"x1"], rect.xxyy[,"x2"])-.5, usr[1]), pmax(pmin(rect.xxyy[,"y1"], rect.xxyy[,"y2"])-.5, usr[3]), pmin(pmax(rect.xxyy[,"x1"], rect.xxyy[,"x2"])+.5, usr[2]), pmin(pmax(rect.xxyy[,"y1"], rect.xxyy[,"y2"])+.5, usr[4]), col=rect.col, border=NA) # rectangle # Reference variables on their background: if(ref.ny>0) rect(blot.xlim[1], blot.ylim[1], usr[2], usr[3], col=ref.col, border=NA) if(ref.nx>0) rect(blot.xlim[1], blot.ylim[1], usr[1], usr[4], col=ref.col, border=NA) # Mark-up strips: for(y in strips.y) if(y>usr[3]+.5) rect(blot.xlim[1], y-.5, blot.xlim[2], y+.5, col=strips.col, border=NA) # hor strips for(x in strips.x) if(x>usr[1]+.5) rect(x-.5, blot.ylim[1], x+.5, blot.ylim[2], col=strips.col, border=NA) # ver strips # Grid for reference variables: if(ref.flag) { ref.x <- blot.xlim[1] + seq(length=ref.nx) - .5 ref.y <- blot.ylim[1] + seq(length=ref.ny) - .5 # x fast, y slow: tab[x,y] if(max(ref.nx,ref.ny)>0) { xrefgrid <- c(rep(ref.x, ref.ny), # bottom left: ref vars rep(nxsel, ref.ny), # bottom: hor ref vars rep(ref.x, length(nysel)) ) # left: ver ref vars yrefgrid <- c(rep(ref.y, rep(ref.nx,ref.ny)), # bottom left: ref vars rep(ref.y, rep(length(nxsel),ref.ny)), # bottom: hor ref vars rep(nysel, rep(ref.nx,length(nysel))) ) # left: ver ref vars reftab <- c(c(tab[ref.vars.x, ref.vars.y]), c(tab[nxsel, ref.vars.y]), c(tab[ref.vars.x, nysel]) ) xrefgrid <- xrefgrid[!is.na(reftab)] yrefgrid <- yrefgrid[!is.na(reftab)] reftab <- reftab[!is.na(reftab)] # Error message 2010/09/19: ## Error in rgb(pmin(1 - reftab, 1), 1 - abs(reftab), pmin(1 + reftab, 1), : ## color intensity nan, not in [0,1] if(color.scale) { # heatmap points(xrefgrid, yrefgrid, col=rgb(pmin(1-reftab,1),1-abs(reftab),pmin(1+reftab,1),maxC=1), pch=blot.pch, cex=blot.sz/max(diff(blot.xlim),diff(blot.ylim))) } else { # blockplot ("blot") refcex <- abs(reftab)/(max(diff(blot.xlim),diff(blot.ylim))/blot.sz) refcol <- if(tab.abs) { blot.col.pos } else { c(blot.col.neg,0,blot.col.pos)[2+sign(reftab)] } points(x=xrefgrid, y=yrefgrid, col=refcol, pch=blot.pch, cex=refcex) } } } # ----------------------------- # Main action: Draw blocks, but draw them as "." if there are < 2pixels per block pixels.per.block.x <- diff(grconvertX(par()$usr[1:2], from="user", to="device")) / # Hor pixel number of plotting area (right-left) diff(par()$usr[1:2]) # Hor number of blocks pixels.per.block.y <- diff(grconvertY(par()$usr[4:3], from="user", to="device")) / # Vertical pixel number of plotting area (bottom-top) diff(par()$usr[3:4]) # Ver number of blocks if(color.scale) { points(xgrid, ygrid, col=rgb(pmin(1-tb,1),1-abs(tb),pmin(1+tb,1),maxC=1), pch=blot.pch, cex=blot.sz/max(diff(blot.xlim),diff(blot.ylim))) } else if(min(pixels.per.block.x, pixels.per.block.y) < 2 & !hdcopy) { points(xgrid, ygrid, col=xycol, pch=".", cex=1) # Draw points!!!!! Main action -- helicopter view, pch="." } else { points(xgrid, ygrid, col=xycol, pch=blot.pch, cex=xycex) # Draw blocks!!!!! Main action -- detailed view, pch=square } # ----------------------------- # Diagonal blocks: Not drawn for #Missing or #Complete Pairs: ## if(length(grep("^[#]",tab.has))==0 & length(nxsel)>0 & length(nysel)>0) { ## diagmin <- max(min(nxsel),min(nysel)) ## diagmax <- min(max(nxsel),max(nysel)) ## if(abs(diagmin)==Inf | abs(diagmax)==Inf) { print(nxsel); print(nysel); cat("-----------\n") } # diagnostic ## if(diagmin0) { lab.bot.left <- bquote(bold(.(paste(tab.has.parts[1], "\n<", pvals.thresh.lst[pvals.thresh.i], "\n", tab.has.parts[2], sep="")))) } else { lab.bot.left <- # bquote(bold(.(paste(tab.has,"\n>=",cors.thresh.lst[ cors.thresh.i],sep="")))) bquote(bold(.(paste(tab.has.parts[1], if(as.numeric(cors.thresh.lst[ cors.thresh.i]) > 0) { paste("\n| |>=", cors.thresh.lst[ cors.thresh.i]) }, "\n", tab.has.parts[2], sep="")))) } if(hdcopy) options(warn=-1) # Ignore warnings; PDF has a problem with bold character 0xa = "\n" probably. ## mtext(side=2, at=1.1*par()$usr[3]-0.1*par()$usr[4], text=lab.bot.left, cex=.8, las=1, adj=1.1) mtext(side=2, at=1.1*par()$usr[3]-0.1*par()$usr[4], text=lab.bot.left, cex=.8, las=1, adj=1.1) if(hdcopy) options(warn=1) # Back to reporting warnings at the end; +1=immediate reporting; 0=reporting after exit } # ------------------------ End 'a.nav.blockplot()' # a.nav.lens.add <- function() { # Add two vars to exam list for lens plot and scatterplot brushing. if(verbose) cat("\na.nav.lens.add") ref.nx <- if(ref.flag) length(ref.vars.x) else 0 # Reference variables on left and bottom ref.ny <- if(ref.flag) length(ref.vars.y) else 0 # Reference variables on left and bottom usr <- c(blot.xlim[1]+ref.nx, blot.xlim[2], blot.ylim[1]+ref.ny, blot.ylim[2]) # Window for non-reference variables if(xhair.xy[1]>=usr[1]) { ix <- round(xhair.xy[1]) } else if(xhair.xy[1]>=blot.xlim[1]) { ix <- ref.vars.x[round(xhair.xy[1]-blot.xlim[1]+.5)] } else { ix <- -1 } if(xhair.xy[2]>=usr[3]) { iy <- round(xhair.xy[2]) } else if(xhair.xy[2]>=blot.ylim[1]) { iy <- ref.vars.y[round(xhair.xy[2]-blot.ylim[1]+.5)] } else { iy <- -1 } if(1<=ix & ix<=ncol(dat) & 1<=iy & iy<=ncol(dat) & mode=="Blockplot") { lens.i <<- length(lens.lst.xy) + 1 lens.lst.xy[[lens.i]] <<- c(x=colnames(dat)[ix],y=colnames(dat)[iy]) } if(length(lens.lst.xy)>=2) # Remove last if same as preceding. if(all(lens.lst.xy[[lens.i]]==lens.lst.xy[[lens.i-1]])) { lens.lst.xy <<- lens.lst.xy[-length(lens.lst.xy)] lens.i <<- lens.i - 1 } a.nav.scatterplot() } # ------------------------ End 'a.nav.lens.add()' # a.nav.splom <- function(equalaxes=F, jitter=.2, extra=.1, ...) { if(verbose) cat("\na.nav.splom") # Find variables sel.vars <- a.nav.su(c(strips.x, strips.y)) if(length(sel.vars)<2) return() x <- rbind(dat[,sel.vars],NA,NA) # Add rows for extremes to stretch the range nxrow <- nrow(x) colnames(x) <- gsub(bl.sep,"\n",colnames(x)) for(j in 1:ncol(x)) { vals <- a.nav.su(x[,j]) nvals <- length(vals) if(nvals<=lens.max.ncat) { space <- min(diff(vals)) x[,j] <- x[,j] + runif(nxrow, -space*jitter, space*jitter) } rg <- range(x[,j], na.rm=T) x[nxrow-c(1,0),j] <- mean(rg) + c(-1,1)*diff(rg)/2*(1+extra) } if(equalaxes) x[nxrow-c(1,0),] <- range(x, na.rm=T) # Cyclic repetition a.nav.to.splom.window() par(mgp=blot.mgp, mar=blot.mar) pairs(x, col=c(col.case,NA,NA), gap=splom.gap, pch=splom.pch, row1attop=F, cex=splom.cex, cex.labels=splom.cex.labels, cex.axis=splom.cex.axis, ...) a.nav.to.main.window() } # ------------------------ End 'a.nav.splom()' # a.nav.scatterplot <- function(hdcopy=F, jitter=.2, extra=.1, cex.lab=1, cex.axis=.8, mgp=c(1.8,.5,0), mar=c(4.5,4.5,1,1), cex=0.8, pch=16, ...) { if(verbose) cat("\na.nav.scatterplot") # Sanity checks: if(length(lens.lst.xy)==0) { cat("--- a.nav.scatterplot: lens.lst.xy not set\n"); return() } if(length(lens.i)==0) { cat("--- a.nav.scatterplot: lens.i not set\n"); return() } if(lens.i<1 | lens.i>length(lens.lst.xy)) { cat("--- scatter: lens.i out of bound:",lens.i,"\n"); return() } if(!all(lens.lst.xy[[lens.i]] %in% colnames(dat))) { cat("--- scatter: lens variables not in colnames(dat):",lens.lst.xy[[lens.i]],"\n"); return() } if(lens.hor.vert) { xyvar <- lens.lst.xy[[lens.i]] } else { xyvar <- rev(lens.lst.xy[[lens.i]]) } # Following is for jitter, axes, traces/curves: xy <- dat[,xyvar] # Keeps the column labels if(adjust.flag) for(j in 1:2) if(xyvar[j] %in% colnames(adjust.cors)) xy[,j] <- adjust.dat[,xyvar[j]] # Pre-jitter, complete axes, not just for complete pairs: xlim <- a.nav.stretch(range(xy[,1], na.rm=T), 1.05) ylim <- a.nav.stretch(range(xy[,2], na.rm=T), 1.05) # Needed to refer back to dat: sel.xy <- !is.na(xy[,1]) & !is.na(xy[,2]) if(sum(sel.xy)<=1) return() # Continue only if there are more than two complete pairs: xy.sub <- xy.sub.raw <- xy[sel.xy,] # Raw for calculating correlation sux <- a.nav.su(xy.sub[,1]) suy <- a.nav.su(xy.sub[,2]) nxvals <- length(sux) nyvals <- length(suy) # Mean trace/smooth curve: Needs to be computed before jittering. if(curve.draw) { if(nxvals<=lens.max.ncat) { tmp <- tapply(xy.sub[,2], xy.sub[,1], FUN=mean) decor.trace <- cbind(x=as.numeric(names(tmp)), y=c(tmp)) } else { decor.trace <- lowess(xy.sub) } } # Labels: xlab <- bquote(bold(.(colnames(xy.sub)[1]))) ylab <- bquote(bold(.(colnames(xy.sub)[2]))) # Axes: depending on number of values xaxt <- ifelse(nxvals<=lens.max.ncat, "n", "s") yaxt <- ifelse(nyvals<=lens.max.ncat, "n", "s") # Jitter the two variables if fewer than 'lens.max.ncat' values if(nxvals<=lens.max.ncat & length(sux)>1) { space <- min(diff(sux),na.rm=T) xy.sub[,1] <- xy.sub[,1] + runif(nrow(xy.sub), -space*jitter, space*jitter) xlim <- xlim + c(-1,1)*space*jitter } if(nyvals<=lens.max.ncat & length(suy)>1) { space <- min(diff(suy),na.rm=T) xy.sub[,2] <- xy.sub[,2] + runif(nrow(xy.sub), -space*jitter, space*jitter) ylim <- ylim + c(-1,1)*space*jitter } # # Plotting in main window: if(mode=="Scatterplot") { if(!hdcopy) a.nav.to.main.window() par(mar=mar, mgp=mgp) # mar, mgp from a.nav.scatterplot() defaults, not from the av environment # Brushing points if on device 1: Needs to be done before jittering if(any(is.na(brush))) { usr <- par()$usr brush <<- c(mean(usr[1:2]), mean(usr[1:2])+diff(usr[1:2])/10, mean(usr[3:4]), mean(usr[3:4])+diff(usr[3:4])/10) } if(mouse.button==0) { # Move-brushing brush <<- c(mouse.xy[1], mouse.xy[1] + brush[2] - brush[1], mouse.xy[2], mouse.xy[2] + brush[4] - brush[3]) } else if(mouse.button==2) { # Resize-brushing brush[c(1,3)] <<- c(mouse.xy[1], mouse.xy[2]) }## else ## if(mouse.button==2 & brush.reset=="Reset") { # Reset brush ## usr <- par()$usr ## brush <<- c(mouse.xy[1], mouse.xy[1]+.02*(usr[2]-usr[1]), mouse.xy[2], mouse.xy[2]+.02*(usr[4]-usr[3])) ## brush.reset <<- "" ## } # Set colors under brush: sel <- which((min(brush[1:2]) <= xy.sub[,1]) & (xy.sub[,1] <= max(brush[1:2])) & (min(brush[3:4]) <= xy.sub[,2]) & (xy.sub[,2] <= max(brush[3:4])) ) if(brush.col!="gray") col.case[sel.xy][sel] <<- brush.col # "gray" = non-brushing color, for counting only plot(xy.sub, xlim=xlim, ylim=ylim, xlab=xlab, ylab=ylab, cex=cex, pch=pch, col=col.case[sel.xy], xaxt=xaxt, yaxt=yaxt, ...) if(xaxt=="n") { axis(side=1, at=sux) } if(yaxt=="n") { axis(side=2, at=suy) } # Show brush if on device 1 and count under brush; else show statistics and histograms rect(brush[1], brush[3], brush[2], brush[4], border=brush.col, col=NULL) if(!hdcopy) { mtext(side=1, line=1, at=1.01*par()$usr[1]-0.01*par()$usr[2], text="#Brushed:", cex=.6, adj=1) mtext(side=1, line=2, at=1.01*par()$usr[1]-0.01*par()$usr[2], text=length(sel), cex=.6, adj=1) } } # End if(mode=="Scatterplot") # Plotting in lens window: a.nav.to.lens.window() par(mfrow=c(3,1), mar=mar, mgp=mgp) plot(xy.sub, xlim=xlim, ylim=ylim, xlab=xlab, ylab=ylab, cex=cex, pch=pch, col=col.case[sel.xy], xaxt=xaxt, yaxt=yaxt, ...) if(xaxt=="n") { axis(side=1, at=sux) } if(yaxt=="n") { axis(side=2, at=suy) } # Smooth or mean trace: if(curve.draw) { lines(decor.trace, lwd=curve.lwd, col=curve.col) if(nxvals<=lens.max.ncat) points(decor.trace, col=curve.col, pch=16, cex=curve.cex) } # Report correlation, p-value, N: mtext(text=bquote(bold(paste("Frame #",.(lens.i),": n = ",.(nrow(xy.sub)),sep=""))), side=1, line=4, cex=.8) pval <- round(2*pnorm(abs(cor(xy.sub.raw[,1],xy.sub.raw[,2])),m=0,s=1/sqrt(notna[xyvar[1],xyvar[2]]),low=F),8) ## switch(pval.type, ## "perm" = round((pvals[ix,iy]+pvals.nsim)/2/pvals.nsim,8), ## "standard" = round(2*pnorm(abs(cors.compl[ix,iy]),m=0,s=1/sqrt(notna[ix,iy]),low=F),8) ## ) mtext(text=bquote(bold(paste("Corr = ",.(round(cor(xy.sub.raw[,1],xy.sub.raw[,2]),3)), " (pval = ",.(pval),")",sep=""))), #xxxxxxxx # 2*pnorm(abs(cors.imp),m=0,s=1/sqrt(notna),low=F) # standard p-values side=1, line=5.5, cex=.8) # Plot histograms: par(mar=mar+c(0,0,3,0), mgp=mgp) # Space above for statistics xlab <- bquote(bold(.(xyvar[1]))) hist(xy[,1], xlim=xlim, xlab=xlab, xaxt=xaxt, breaks=20, col="gray", cex.lab=1, cex.axis=.8, main="") if(xaxt=="n") { axis(side=1, at=sux) } par(mar=mar+c(3,0,0,0), mgp=mgp) ylab <- bquote(bold(.(xyvar[2]))) # Space below making up for space above in 1st hist hist(xy[,2], xlim=ylim, xlab=ylab, xaxt=yaxt, breaks=20, col="gray", cex.lab=1, cex.axis=.8, main="") if(yaxt=="n") { axis(side=1, at=suy) } a.nav.to.main.window() # To empty plotting queue; else the plot is drawn with delay. } # -------------------------------- End 'a.nav.scatterplot()' # a.nav.xhair.draw <- function(col="gray30") { if(verbose) cat("\na.nav.xhair.draw") lines(x=c(rep(xhair.xy[1],2), NA, xhair.xy[1]+c(-1,1)*diff(blot.xlim)*xhair.sz), y=c(xhair.xy[2]+c(-1,1)*diff(blot.ylim)*xhair.sz, NA, rep(xhair.xy[2],2)), col=col) } # ------------------------ End 'a.nav.xhair.draw()' # a.nav.xhair.move <- function(xy, col="gray30") { if(verbose) cat("\na.nav.xhair.move") a.nav.xhair.draw(col="white") xhair.xy <<- xy a.nav.xhair.draw() } # ------------------------ End 'a.nav.xhair.move()' a.nav.dev.new <- function(...) { ## cat("\na.nav.dev.new: (entry) dev.list() =", dev.list(), "\n") ## cat("\na.nav.dev.new: (entry) dev.list() =", dev.list(), "\n") ## cat("\na.nav.dev.new: platform =", .Platform$OS.type,"\n") if(.Platform$OS.type=="windows") { windows(...) ## cat("\na.nav.dev.new: should have created a new 'windows()'\n") } else if(.Platform$OS.type=="unix") { X11(..., type="Xlib") ## cat("\na.nav.dev.new: should have created a new 'X11()'\n") } ## cat("\na.nav.dev.new: (exit) dev.list() =", dev.list(), "\n") } # ------------------------ End 'a.nav.dev.new()' a.nav.to.main.window <- function() { if(!(dev["main"] %in% dev.list()) ) { a.nav.dev.new(width=15, height=15, title="BLOCKPLOT ") dev["main"] <<- dev.cur() } if(dev.cur()!=dev["main"]) dev.set(dev["main"]) # Minimize dev.set() calls! They make the plot flicker! par(mfrow=c(1,1), mar=blot.mar) ## cat("\na.nav.to.main.window:",dev,"|",dev.list(),"\n") } # ------------------------ End 'a.nav.to.main.window()' a.nav.to.lens.window <- function() { if(!(dev["lens"] %in% dev.list()) ) { # | (dev["lens"] %in% dev[-2]) ) { a.nav.dev.new(width=3.5, height=10, title="LENS ") dev["lens"] <<- dev.cur() } if(dev.cur()!=dev["lens"]) dev.set(dev["lens"]) # cat("\na.nav.to.lens.window:",dev,"|",dev.list(),"\n") } # ------------------------ End 'a.nav.to.lens.window()' a.nav.to.splom.window <- function() { if(!(dev["splom"] %in% dev.list()) ) { # | (dev["splom"] %in% dev[-3]) ) { a.nav.dev.new(width=8, height=8, title="SCATTERPLOT MATRIX ") dev["splom"] <<- dev.cur() } if(dev.cur()!=dev["splom"]) dev.set(dev["splom"]) ## cat("\na.nav.to.splom.window:",dev,"|",dev.list(),"\n") } # ------------------------ End 'a.nav.to.splom.window()' a.nav.in.rects <- function(xy, rects, marg=.5) { # Returns integer() if xy pair is in none of the rects or len(rects)==0 if(is.null(rects) | length(xy)!=2 | length(rects)==0) return(integer()) if(is.null(dim(rects))) rects <- matrix(rects,nrow=1) # It's a vector: one rectangle xx <- cbind(pmin(rects[,1], rects[,2])-marg, pmax(rects[,1], rects[,2])+marg) yy <- cbind(pmin(rects[,3], rects[,4])-marg, pmax(rects[,3], rects[,4])+marg) which(xx[,1]<=xy[1] & xy[1]<=xx[,2] & yy[,1]<=xy[2] & xy[2]<=yy[,2]) } a.nav.stretch <- function(x, fac) { mn <- mean(x); mn + (x-mn)*fac } a.nav.su <- function(x) sort(unique(x)) #================ # Assorted data for colors and help: a.nav.pastell.cols <- c( "aliceblue", "antiquewhite", "antiquewhite1", "antiquewhite2", "azure1", "azure2", "beige", "bisque", "bisque1", "blanchedalmond", "cornsilk", "cornsilk1", "cornsilk2", "darkseagreen1", "darkslategray1", "floralwhite", "gainsboro", "ghostwhite", "gray90", "gray91", "gray92", "gray93", "gray94", "gray95", "gray96", "gray97", "gray98", "honeydew", "honeydew1", "honeydew2", "ivory2", "khaki1", "lavender", "lavenderblush", "lavenderblush1", "lemonchiffon", "lemonchiffon1", "lemonchiffon2", "lightcyan", "lightcyan1", "lightgoldenrod1", "lightgoldenrodyellow", "lightyellow", "lightyellow1", "lightyellow2", "linen", "mintcream", "mistyrose", "mistyrose1", "moccasin", "navajowhite", "navajowhite1", "oldlace", "palegoldenrod", "palegreen", "palegreen1", "paleturquoise", "paleturquoise1", "papayawhip", "peachpuff", "peachpuff1", "seashell", "seashell1", "seashell2", "snow", "snow1", "snow2", "thistle1", "wheat", "wheat1", "wheat2", "whitesmoke") # 1) presentations to SFARI: graphical as much as possible, problem of PDFs # explain "blockplots" as alternatives to heat maps, color coding ~ finance # 2) big hammer: AN # 3) show helicopter view: quit, show size of data in R, restart # 4) zoom around, explain diagonal blocks, instruments # 5) zoom in to ADOS-3: # explain labels, prefixes # point out theory behind labels # highlight rectangle on a_ variables # show how communications variables (high-a_) don't relate to language (low-a_) # instead communications relates to social variables # 6) zoom out some, note the blank correlations ==> missings # show missing value plot # 7) strip search with ADOS PCA # 8) p-values: color coding, 1-pvalue, thresholding (should also provide (-log(1-pval)) # 9) # # # # # a.nav.help <- c(#"________________________________________________________________", " ", " ZOOM AND PAN:", " L-mouse click: Place cross-hair at mouse.", " L-mouse drag: Drag the blockplot ('panning').", " Hit i,I,o,O,=,+,-,_: Zoom in/out wrt cross-hair (caps: faster).", " Hit u,U: Unequal zoom in/out wrt cross-hair; u: xy=-+, U: xy=+-", " Hit arrows: Pan one block left/right/up/down.", " Hit D,d: Pan up/down diagonally.", " Hit space: Accelerated pan; five blocks in direction of last arrow/D/d.", " Hit .,{,},[,]: Pan cross-hair to center/top,bottom/left,right.", " Hit H,V,@: String search a x/y/xy variable to mark with strip and pan to center", " Hit W,w: Forward/backward in window pan/zoom history", "________________________________________________________________", " ", " AESTHETIC ELEMENTS:", " Hit f,F: Font size -,+", " Hit b,B: Block size -,+", " Hit m,M: Margin size; xy=--,++ (same changes)", " Hit n,N: Margin size; xy=-+,+- (opposite changes)", " Hit c,C: Cross-hair size -,+", " Hit t: Toggle blockplot type between color scale and blocksize", "________________________________________________________________", " ", " LENSE WITH SCATTERPLOT AND HISTOGRAMS:", " Hit x: Show variable pair at cross-hair in Scatterplot+Histograms.", " Hit y: Switch x and y in Scatterplot.", " Hit l: Add smooth curve or trace of group means to Scatterplot.", " Hit PgUp,PgDn,Home,End,Del: Lens history - step/jump/delete.", "________________________________________________________________", " ", " SCATTERPLOT BRUSHING IN MAIN WINDOW:", " Hit s: Toggle between Blockplot and Scatterplot in main window.", " In Scatterplot, color the points with a rectangle brush as follows:", " L-mouse drag: Move brush.", " R-mouse drag: Reshape brush.", " Hit S: Cycle colors (scatterplot).", " Hit r: Reset brush size (scatterplot).", "________________________________________________________________", " ", " MARK-UP STRIPS, HORIZONTAL AND VERTICAL:", " R-mouse click/drag: Mark variable(s) with strip(s).", " Hit H,V,@: String search a x/y/xy variable name; mark with strip and pan to center", " Hit z: Show variables marked with strips in scatterplot matrix.", " Hit Z: Same, but with equal axes.", " Hit R: Turn variables marked with strips into reference variables.", " Hit r: Toggle display of reference variables (blockplot).", "________________________________________________________________", " ", " HIGHLIGHT RECTANGLES:", " Hit 1: Highlight rectangle - set corner 1 at cross-hair.", " Hit 2: Highlight rectangle - set corner 2 at cross-hair.", " Hit 3: Highlight rectangle = variable blocks at cross-hair.", " Hit 4: Highlight rectangle = fit rectangle to plotting area (keep aspect ratio).", " Hit 5: Highlight rectangle = fit rectangle to plotting area (ignore aspect ratio).", " Hit 0: Remove highlight rectangles under cross-hair.", "________________________________________________________________", " ", " BLOCKPLOT -- WHAT IT SHOWS:", " Hit ctrl-O: Correlations (complete pairs)", " Hit ctrl-I: Inference (1 - P-values, normal, signed, thresholded)", " Hit <,>: Change threshold on correlations and P-values ('>' for more stringency)", " Hit ctrl-T: Trim extremes from both sides of all variables", " Hit ctrl-A: Toggle between absolute values and signed values", " Hit ctrl-R: Random correlations (all columns randomized)", " Hit ctrl-B: Random correlations (intra-block preserved)", " Hit ctrl-N,ctrl-M: Number of complete/missing pairs for correlation", "________________________________________________________________", " ", " Hit A: Adjustment menu and dialog", "________________________________________________________________", " ", " PRINTING:", " Hit P: Print current plot to file 'blockplotnnn.pdf'.", " Hit J: Print current plot to file 'blockplotnnn.jpg'.", " Hit G: Print current plot to file 'blockplotnnn.png'.", "________________________________________________________________", " ", " RESETS, HELP AND QUIT:", " Hit ctrl-X: Refit Blockplot to window.", " Hit ctrl-C: Clean Blockplot of rectangles.", " Hit ctrl-D: Clean Blockplot of strips.", " Hit h: Help -- What you are looking at.", " Hit ?: Info about data.", " Hit Q (capital!): Quit and return to the R interpreter.", "________________________________________________________________" ) # a.nav.info <- NULL a.nav.info <- c("________________________________________________________________", "", "abc ABC (Aberrant Behavior Checklist)", "abcl-18-59 ABCL 18 - 59 (Adult Behavior Checklist)", "adi-r ADI-R (Autism Diagnostic Interview - Revised)", "adm ADM (Autism Dysmorphology Measure Work Sheet)", "adm2 ADM-2 (Autism Dysmorphology Measure Work Sheet)", "ados-1 ADOS - Module 1 (Autism Diagnostic Observation Schedule)", "ados-2 ADOS - Module 2 (Autism Diagnostic Observation Schedule)", "ados-3 ADOS - Module 3 (Autism Diagnostic Observation Schedule)", "ados-4 ADOS - Module 4 (Autism Diagnostic Observation Schedule)", "bapq BAPQ (The Broad Autism Phenotype Questionnaire)", "cbcl-2-5 CBCL 2 - 5 (Child Behavior Checklist)", "cbcl-6-18 CBCL 6 - 18 (Child Behavior Checklist)", "ctopp CTOPP (Comprehensive Test of Phonological Processing)", "ctrf-2-5 C-TRF (Caregiver / Teacher Report Form for ages 1.5 to 5 years)", "das-ii-early-years DAS-II - Early Years (Differential Ability Scales, 2nd Edition)", "das-ii-school-age DAS-II - School Age (Differential Ability Scales, 2nd Edition)", "dcdq DCDQ (Developmental Coordination Disorder Questionnaire)", "fhi-informant FHI-Informant (Family History Interview)", "fhi-interviewer FHI-I (IMGSAC Interviewer''s Impressions)", "fhi-subject FHI - Subject (Family History Interview)", "mullen Mullen Scales of Early Learning ", "ppvt-4a PPVT-4, Form A (Peabody Picture Vocabulary Test, 4th Edition)", "ppvt-4b PPVT-4, Form B (Peabody Picture Vocabulary Test, 4th Edition)", "purdue-pegboard Purdue Pegboard ", "ravens-coloured Raven''s Coloured Progressive Matrices ", "ravens-standard Raven''s Standard Progressive Matrices ", "rbs-r RBS - R (Repetitive Behavior Scale - Revised)", "srs-adult SRS-ARV (Social Responsiveness Scale - Adult Research Version)", "srs-parent Parent SRS (Social Responsiveness Scale)", "srs-teacher Teacher SRS (Social Responsiveness Scale)", "ssc-background-hx SSC Background History Form ", "ssc-diagnosis SSC Diagnosis Form ", "ssc-hwhc SSC Height, Weight, Head Circumference ", "ssc-med-hx SSC Medical History Form ", "ssc-med-hx-v2 SSC Medical History Form - Version 2 ", "ssc-pedigree SSC Pedigree ", "ssc-psuh SSC Parent Substance Use History ", "ssc-treatment-hx SSC Treatment History Form ", "trf-6-18 TRF 6 - 18 (Teacher Report Form for ages 6 to 18 years)", "vineland-ii Vineland-II (Vineland Adaptive Behavior Scales, 2nd Ed.", "wasi WASI (Wechsler Abbreviated Scale of Intelligence)", "wisc-iv WISC-IV (Wechsler Intelligence Scale for Children, 4th Edition)", "", "________________________________________________________________" ) ## #xxxx type in SRS questionnaire items, the p1_pSRS version; srs.questionnaire, ## a.nav.HelpSRS <- c("________________________________________________________________", ## "", ## "________________________________________________________________" ## ) ## a.nav.HelpADIcoding <- ## c("________________________________________________________________", ## "", ## "0: Behavior of type specified in the coding is or was not present.", ## "1: Behavior is or was present in an abnormal form (or lack of behavior), but not severe or frequent for 2.", ## "2: Definite abnormality meets or met criteria.", ## "3: Definite abnormality, and more severe than 2.", ## "7: Definite abnormality in the general area of the coding, but not of the type specified.", ## "8: Not applicable (because outside the relevant age range, lacks language, or has not shown loss of skills.", ## "9: Not known or not asked.", ## "", ## "993: Milestone achieved, then relapsed", ## "994: Milestone never achieved", ## "995: Milestone still not reached (less than a year)", ## "996: Not known; apparently normal", ## "997: Not known; apparently delayed", ## "998: Not applicable (e.g., physical handicap)", ## "999: Not known or not asked", ## "________________________________________________________________") #================================================================ ## # ASSOCIATION NAVIGATOR: ## # Load this file with: source("association-navigator.R") ## # Load data: generated with ## # with(cv, save(dat, pvals, pvals.nsim, psums, psqrs, file="aut-new.RData")) ## with(cv, qqnorm(psums[col(psums)>row(psums)],pch=".")) ## with(cv, pnorms <- 2*pnorm(abs(cors.imp), ## m=0, #psums/pvals.nsim, ## s=sqrt(psqrs/pvals.nsim), # sdev=sqrt(E(X^2)-E(X)^2) does not work ## lower.tail=F)) ## with(cv, sum(is.na(pnorms))) # [1] 928 ## with(cv, sum(is.na(diag(pnorms)))) # [1] 0: NAs are off=diagonal ## with(cv, sum(is.na(pvals))) # [1] 0 ## with(cv, sum(is.na(cors.imp))) # [1] 0 ## with(cv, sum(cors.imp==0)) # [1] 200 ## with(cv, all(diag(pnorms)==0)) # [1] FALSE ## with(cv, hist(pnorms[col(pnorms)>row(pnorms)], col="gray", breaks=1000)) ## with(cv, {sel <- col(pnorms)>row(pnorms) ## plot(x=(pvals[sel]+pvals.nsim)/2/pvals.nsim, ## y=pnorms[sel], pch=".")} ) ## with(cv, {sel <- col(pnorms)>row(pnorms); plot(c(psums[sel]/pvals.nsim), pch=".")}) ## with(cv, {sel <- col(psqrs)>row(pnorms); plot(c(psqrs[sel]/pvals.nsim), pch=".")}) ## ## load("aut-pvals.RData"); load("aut.RData") ## # Create association viewer for 'aut.dat': ## bv <- a.nav.create(aut.dat, null.grps=aut.dat[,"module_ADOS-"], grp.var=aut.dat[,"release_ID-"]) ## with(bv, { rm(pvals, pvals.nsim); pvals <- pvals; pvals.nsim <- pvals.nsim }) ## # Accumulate p-value simulations as needed: ## ## with(bv, corr.pvals.marg(nsim=100)) ## # This line can be called repeatedly with increasing 'nsim'; it keeps accumulating. ## # Restore saved pvals simulation data: ## ## load("../pvals-et-al.RData") ## ## with(bv, { pvals <- pvals.back ## ## psums <- psums.back ## ## psqrs <- psqrs.back ## ## pvals.nsim <- pvals.nsim.back} ) ## # Run association viewer: ## a.nav.run(bv) ## # Lifted from between functions above: ## with(dv, { ## ix <- round(mouse[1]) ## ixlab <- rownames(tab)[ix] ## ixlab <- rev(strsplit(ixlab,bl.sep)[[1]])[1] ## ixlab <- substring(ixlab,1,nchar(ixlab)-1) ## ixmatches <- grep(ixlab,rownames(tab)) ## ixsets <- cumsum(c(1,diff(ixmatches)!=1)) ## rbind(ixmatches,ixsets) ## range(ixmatches[ixsets==ixsets[ixmatches==ix]], na.rm=T) ## } ) #---------------------------------------------------------------- # Backing up: # cp -R -u -v -p c:/D/* d:/D/ # cp -R -u -v -p "c:/D/D-AUTISM/DATA-2009-08-18-SFARI Collection V5"/* "d:/D/D-AUTISM/DATA-2009-08-18-SFARI Collection V5"/ #---------------------------------------------------------------- # # # To do: # # - Dissemination: # . Send .RData to Brown U # . Write paper/documentation -- done # . Create package # . Find public data example # # - Functionality: # . Diagonal blocks: allow diagonal entries !=1, as in cell count plot; choose slightly different blue/gray color # . Sensitivity analysis: remove 1,2,... extreme obs on both sides and recalculate corrs AND p-vals # (ideas: might estimate sd's (denominator) indep of cov (numerator); remove extremes only in cov, leave sd's alone...?) # . Measures of sensitivity/specificity for hor-ver quadrant divisions in scatterplot # These would be alternatives to correlation. # . Brushing in lens with its own getGraphicsEvent() running? # . Use points() instead of rect() when zoomed out # . Ordering/subselecting variables: reorder whole instruments and individual variables # Most frequently needed is for example to bring along Age, ADOS Module,... # I.e., select variables in order and insert a copy anywhere (enlarging the matrix) # . Removing subsets (e.g., outlier>70 for age_q05a-walked_unaideda_ADI) # . Subgrouping for exclusion, according to variables, generalizing ctrl-G # . Graphical choice menu for colors and glyphs for more choices and avoiding color cycling # . Numeric table representation instead of blocks # . Add numeric values of strongest correlation in the visible window to bottom left labels # . Cross-marking: clicking a label, highlight both row and col (to find corrs of 2 rows) # . Margins 3 and 4 to move hor/ver on 'mousedown' # . Drill-down to individuals and display them (value and extremity for each variable?) # - Statistics: # . Adjusted correlations, e.g., for age or ADOS level, or instrument adjusted for PC1 # for age, language level, ... (user chosen) # . Robust correlations: winsorize? Pull in 'k' most extreme observations on either side. # . Zero-centered correlations # . Permutations of non-NAs: should be doable w/o exorbitant cost # . New methodology: penalized SVD with squared lasso for clustering of cases and vars # - Data: # . Other well-populated tables: dcdq; dcdq.raw; ppvt.4a # . Medical v2 tables # . Master table with corrs, pvals --> Daughter tables with their corrs, pvals # . General arrangement of variables: Summaries - PCA - items (check ADOS, ADI, ...) # . Sort items by PC loadings for instruments with non-organized items. # General concept: Sort by descending correlation with a given set of variables, # but only within instrument? # Global does not work: needs to be within instrument # with(bv, { v <- c("PC1_ADOS2-","PC2_ADOS2-","PC3_ADOS3-") # ord <- order(apply(tab[,v], 1, function(x) max(abs(x))), decreasing=T) # tab <- tab[ord,ord] # }) # . New dataset consisting of union of PROBS and SIBS for Vineland, CBCL, pSRS, tSRS # . New dataset for difference prob-sib for Vineland, CBCL, pSRS, tSRS # - New Tools: # . New blockplot navigator for raw data, rectangular, that permits sorting by columns # . Other lenses: parallel coordinate lense # histogram/barplot array of marked variables # # - Done: # . Re-engineered away the use of the middle mouse button from both blockplot and scatterplot 2010/11/14 # . ctrl-M for number of missings in cell; complement of ctrl-N (actually: ctrl-M shows as ctrl-J) # . Find-tool to zoom to variables known approximately by name: # enter "string_hor" and/or "string_ver" by popup, # grep("string",colnames(dat),val=T), # offer choices by menu when non-unique, # then center hor and/or ver accordingly # . Reference variables in the left & lower margins !!!!! # . Need a mode that shows magnitude without sign (ctrl-A) # . SRS: merge labels to 'p1pSRS', move 's1xSRS' to end, so p1xSRS are contiguous # . Dataset of summary variables 'aut.smry' in addition to 'aut.dat' # . Aspect zoom: use u/U (next to i/I and o/O) to shrink/stretch aspect ratio # . Unequal margins: use n/N (next to m/M) to widen one margin and narrow the other # . Vineland outlier in 'composite_std_score_p1' # . Other lenses: scatterplot matrix of marked variables # . Statistics: Normal p-values from permutations: means & sd # . Random correlations based on blocked permutations DONE 2010/04/20 or so # . Lens: if N=0, don't draw, ignore; currently crashes 2010/07 # . PCs: add '% accounted for' to label, make it 4 PCs, not 3, in all tables # . Add more PCs to ABC # 2010/08/12: probably from 'african-american.fa_raceParent' ## Warning messages: ## 1: In min(diff(vals)) : no non-missing arguments to min; returning Inf ## 2: In runif(nxrow, -space * jitter, space * jitter) : NAs produced ## 3: In min(x) : no non-missing arguments to min; returning Inf ## 4: In max(x) : no non-missing arguments to max; returning -Inf ## 5: In min(diff(vals)) : no non-missing arguments to min; returning Inf ## 6: In runif(nxrow, -space * jitter, space * jitter) : NAs produced ## 7: In min(x) : no non-missing arguments to min; returning Inf ## 8: In max(x) : no non-missing arguments to max; returning -Inf # Hit "V", "H", "@": Replace menu(...) with select.list(..., multiple=T) # (Go to first selected item in list if multiple chosen.) # (Use preselect=... from previous if overlap?) # Add srs.questionnaire to info. Try to control the font. # Debug the following which occurs after 'ctrl-T' (blockwise random corrs): ## Error in apply(dat[, postfx == pf], 1, function(x) all(is.na(x))) : ## dim(X) must have a positive length ## In addition: Warning message: ## In cor(dat.sim, use = "pairwise.complete.obs") : ## the standard deviation is zero # Debug this message: after erroneously hitting 4 w/o 3 beforehand ## Warning messages: ## 1: In min(nxsel) : no non-missing arguments to min; returning Inf ## 2: In min(nysel) : no non-missing arguments to min; returning Inf ## 3: In max(nxsel) : no non-missing arguments to max; returning -Inf ## 4: In max(nysel) : no non-missing arguments to max; returning -Inf ## 5: In min(nxsel) : no non-missing arguments to min; returning Inf ## 6: In min(nysel) : no non-missing arguments to min; returning Inf ## 7: In max(nxsel) : no non-missing arguments to max; returning -Inf ## 8: In max(nysel) : no non-missing arguments to max; returning -Inf # Warnings: ## Warning messages: ## 1: In min(nysel) : no non-missing arguments to min; returning Inf ## 2: In max(nysel) : no non-missing arguments to max; returning -Inf ## 3: In min(nysel) : no non-missing arguments to min; returning Inf ## 4: In max(nysel) : no non-missing arguments to max; returning -Inf ## ... #================================================================