Métodos estatísticos para plotar dados de maneira mais eficiente quando milhões de pontos estão presentes?

Acho que R pode levar muito tempo para gerar gráficos quando milhões de pontos estão presentes - o que não é surpreendente, uma vez que os pontos são plotados individualmente. Além disso, essas parcelas costumam ser muito confusas e densas para serem úteis. Muitos dos pontos se sobrepõem e formam uma massa negra e é gasto muito tempo plotando mais pontos nessa massa.

Existem alternativas estatísticas para representar grandes dados em um gráfico de dispersão padrão? Eu considerei um gráfico de densidade, mas que outras alternativas existem?$n$

Essa é uma tarefa difícil, sem soluções prontas (é claro, porque o gráfico de densidade é um substituto tentador do que ninguém realmente se importa). Então o que você pode fazer?

Se eles realmente se sobrepuserem (ou seja, tiverem exatamente as mesmas coordenadas X e Y) e você não estiver usando alfa, a melhor idéia seria apenas reduzir a sobreposição usando unique(com alfa, pode ser somado a esses grupos).

Caso contrário, você pode arredondar manualmente as coordenadas para os pixels mais próximos e usar o método anterior (ainda assim, essa é uma solução suja).

Finalmente, você pode fazer um gráfico de densidade apenas para usá-lo para subamostrar os pontos nas áreas mais densas. Por outro lado, isso não fará exatamente o mesmo gráfico e poderá introduzir artefatos se não for ajustado com precisão.

Veja o pacote hexbin que implementa o papel / método de Dan Carr. A vinheta em pdf tem mais detalhes que cito abaixo:

1. Visão Geral

Binning hexagonal é uma forma de histograma bivariado útil para visualizar a estrutura em conjuntos de dados com n grande. O conceito subjacente de bin hexagonal é extremamente simples;

1. o plano xy sobre o conjunto (intervalo (x), intervalo (y)) é mosaico por uma grade regular de hexágonos.
2. o número de pontos que caem em cada hexágono é contado e armazenado em uma estrutura de dados
3. os hexágonos com contagem> 0 são plotados usando uma rampa de cores ou variando o raio do hexágono na proporção das contagens. O algoritmo subjacente é extremamente rápido e eficaz para exibir a estrutura de conjuntos de dados com$n \ge 10^6$

Se o tamanho da grade e os cortes na rampa de cores forem escolhidos de maneira inteligente, a estrutura inerente aos dados deverá emergir nos gráficos em bin. As mesmas advertências aplicam-se ao escaneamento hexagonal e aos histogramas e deve-se tomar cuidado ao escolher os parâmetros de escaneamento

Devo admitir que não compreendo completamente seu último parágrafo:

"Não estou procurando um gráfico de densidade (embora esses sejam frequentemente úteis), gostaria da mesma saída que uma chamada de gráfico simples, mas muito mais rápido do que milhões de overplots, se possível."

Também não está claro que tipo de gráfico (função) você está procurando.

Dado que você tem variáveis ​​métricas, poderá encontrar gráficos em hexágono ou em girassol úteis. Para mais referências, consulte

• Gráficos de grandes conjuntos de dados by Unwin/Theus/Hofmann
• Quick-R on "High Density Scatterplots"
• ggplot2's stat_hexbin

Another direct answer to the question is the rgl package, which can plot millions of points using OpenGL. Also, specify a point size (e.g. 3) and zoom out to see these centers of masses as monolithic blocks, or zoom in and see the structure of what used to be monolithic - the point sizes are constant but the distances among them on the screen depend on the zooming. Alpha levels can also be used.

Here's a file I call bigplotfix.R. If you source it, it will define a wrapper for plot.xy which "compresses" the plot data when it is very large. The wrapper does nothing if the input is small, but if the input is large then it breaks it into chunks and just plots the maximum and minimum x and y value for each chunk. Sourcing bigplotfix.R also rebinds graphics::plot.xy to point to the wrapper (sourcing multiple times is OK).

Note that plot.xy is the "workhorse" function for the standard plotting methods like plot(), lines(), and points(). Thus you can continue to use these functions in your code with no modification, and your large plots will be automatically compressed.

This is some example output. It's essentially plot(runif(1e5)), with points and lines, and with and without the "compression" implemented here. The "compressed points" plot misses the middle region due to the nature of the compression, but the "compressed lines" plot looks much closer to the uncompressed original. The times are for the png() device; for some reason points are much faster in the png device than in the X11 device, but the speed-ups in X11 are comparable (X11(type="cairo") was slower than X11(type="Xlib") in my experiments).

A razão pela qual escrevi isso é porque estava cansado de correr plot()acidentalmente em um grande conjunto de dados (por exemplo, um arquivo WAV). Nesses casos, eu teria que escolher entre aguardar vários minutos para a plotagem terminar e encerrar minha sessão R com um sinal (perdendo, assim, meu histórico de comandos e variáveis ​​recentes). Agora, se me lembro de carregar esse arquivo antes de cada sessão, posso obter uma plotagem útil nesses casos. Uma pequena mensagem de aviso indica quando os dados da plotagem foram "compactados".

# bigplotfix.R
# 28 Nov 2016

# This file defines a wrapper for plot.xy which checks if the input
# data is longer than a certain maximum limit. If it is, it is
# downsampled before plotting. For 3 million input points, I got
# speed-ups of 10-100x. Note that if you want the output to look the
# same as the "uncompressed" version, you should be drawing lines,
# because the compression involves taking maximum and minimum values
# of blocks of points (try running test_bigplotfix() for a visual
# explanation). Also, no sorting is done on the input points, so
# things could get weird if they are out of order.
test_bigplotfix = function() {
  oldpar=par();
  par(mfrow=c(2,2))
  n=1e5;
  r=runif(n)
  bigplotfix_verbose<<-T
  mytitle=function(t,m) { title(main=sprintf("%s; elapsed=%0.4f s",m,t["elapsed"])) }
  mytime=function(m,e) { t=system.time(e); mytitle(t,m); }

  oldbigplotfix_maxlen = bigplotfix_maxlen
  bigplotfix_maxlen <<- 1e3;

  mytime("Compressed, points",plot(r));
  mytime("Compressed, lines",plot(r,type="l"));
  bigplotfix_maxlen <<- n
  mytime("Uncompressed, points",plot(r));
  mytime("Uncompressed, lines",plot(r,type="l"));
  par(oldpar);
  bigplotfix_maxlen <<- oldbigplotfix_maxlen
  bigplotfix_verbose <<- F
}

bigplotfix_verbose=F

downsample_xy = function(xy, n, xlog=F) {
  msg=if(bigplotfix_verbose) { message } else { function(...) { NULL } }
  msg("Finding range");
  r=range(xy$x);
  msg("Finding breaks");
  if(xlog) {
    breaks=exp(seq(from=log(r[1]),to=log(r[2]),length.out=n))
  } else {
    breaks=seq(from=r[1],to=r[2],length.out=n)
  }
  msg("Calling findInterval");
  ## cuts=cut(xy$x,breaks);
  # findInterval is much faster than cuts!
  cuts = findInterval(xy$x,breaks);
  if(0) {
    msg("In aggregate 1");
    dmax = aggregate(list(x=xy$x, y=xy$y), by=list(cuts=cuts), max)
    dmax$cuts = NULL;
    msg("In aggregate 2");
    dmin = aggregate(list(x=xy$x, y=xy$y), by=list(cuts=cuts), min)
    dmin$cuts = NULL;
  } else { # use data.table for MUCH faster aggregates
    # (see http://stackoverflow.com/questions/7722493/how-does-one-aggregate-and-summarize-data-quickly)
    suppressMessages(library(data.table))
    msg("In data.table");
    dt = data.table(x=xy$x,y=xy$y,cuts=cuts)
    msg("In data.table aggregate 1");
    dmax = dt[,list(x=max(x),y=max(y)),keyby="cuts"]
    dmax$cuts=NULL;
    msg("In data.table aggregate 2");
    dmin = dt[,list(x=min(x),y=min(y)),keyby="cuts"]
    dmin$cuts=NULL;
    #  ans = data_t[,list(A = sum(count), B = mean(count)), by = 'PID,Time,Site']
  }
  msg("In rep, rbind");
  # interleave rows (copied from a SO answer)
  s <- rep(1:n, each = 2) + (0:1) * n
  xy = rbind(dmin,dmax)[s,];
  xy
}

library(graphics);
# make sure we don't create infinite recursion if someone sources
# this file twice
if(!exists("old_plot.xy")) {
  old_plot.xy = graphics::plot.xy
}

bigplotfix_maxlen = 1e4

# formals copied from graphics::plot.xy
my_plot.xy = function(xy, type, pch = par("pch"), lty = par("lty"),
  col = par("col"), bg = NA, cex = 1, lwd = par("lwd"),
  ...) {

  if(bigplotfix_verbose) {
    message("In bigplotfix's plot.xy\n");
  }

  mycall=match.call();
  len=length(xy$x)
  if(len>bigplotfix_maxlen) {
    warning("bigplotfix.R (plot.xy): too many points (",len,"), compressing to ",bigplotfix_maxlen,"\n");
    xy = downsample_xy(xy, bigplotfix_maxlen, xlog=par("xlog"));
    mycall$xy=xy
  }
  mycall[[1]]=as.symbol("old_plot.xy");

  eval(mycall,envir=parent.frame());
}

# new binding solution adapted from Henrik Bengtsson
# https://stat.ethz.ch/pipermail/r-help/2008-August/171217.html
rebindPackageVar = function(pkg, name, new) {
  # assignInNamespace() no longer works here, thanks nannies
  ns=asNamespace(pkg)
  unlockBinding(name,ns)
  assign(name,new,envir=asNamespace(pkg),inherits=F)
  assign(name,new,envir=globalenv())
  lockBinding(name,ns)
}
rebindPackageVar("graphics", "plot.xy", my_plot.xy);

Talvez eu seja evitado pelo meu método, as lembranças ruins de um dos meus profissionais de pesquisa gritando com as pessoas por jogar fora bons dados, traduzindo-os em categorias (é claro, eu concordo hoje em dia lol), não sei. Enfim, se você está falando de um gráfico de dispersão, tive os mesmos problemas. Agora, quando tenho dados numéricos, não faz muito sentido categorizá-los para análise. Mas visualizar é uma história diferente. O que eu achei que funciona melhor para mim é primeiro (1) dividir sua variável independente em grupos usando a função cut. Você pode brincar com o número de grupos e, em seguida, (2) simplesmente plotar o DV contra a versão cortada do IV. R gerará gráficos de caixa em vez daquele gráfico de dispersão nojento. Eu recomendo remover os outliers do plot (use a opção outline = FALSE no comando plot). Novamente, NUNCA desperdiçaria dados numéricos perfeitamente bons categorizando e analisando. Muitos problemas fazendo isso. Embora eu saiba que é um assunto delicado de debate. Mas fazer isso especificamente com o objetivo de, pelo menos, extrair algum sentido visual dos dados, não muito mal que já vi deles. Plotamos dados de até 10 milhões e ainda assim consegui extraí-los com esse método. Espero que ajude! Cumprimentos! Eu vi isso. Plotamos dados de até 10 milhões e ainda assim consegui extraí-los com esse método. Espero que ajude! Cumprimentos! Eu vi isso. Plotamos dados de até 10 milhões e ainda assim consegui extraí-los com esse método. Espero que ajude! Cumprimentos!

For large time series, I have grown to love smoothScatter (part of base R no less). I often have to include some additional data, and preserving the basic plot API is really helpful, for instance:

set.seed(1)
ra <- rnorm(n = 100000, sd = 1, mean = 0)
smoothScatter(ra)
abline(v=25000, col=2)
text(25000, 0, "Event 1", col=2)

Which gives you (if you pardon the design):

It's always available and works well with enormous datasets, so it's nice to at least take a look at what you have.