¿Cómo ajustar una curva suave a mis datos en R?

Estoy tratando de dibujar una curva suave R. Tengo los siguientes datos simples de juguetes:

> x
 [1]  1  2  3  4  5  6  7  8  9 10
> y
 [1]  2  4  6  8  7 12 14 16 18 20

Ahora, cuando lo trazo con un comando estándar, se ve irregular y nervioso, por supuesto:

> plot(x,y, type='l', lwd=2, col='red')

¿Cómo puedo suavizar la curva para que los 3 bordes se redondeen utilizando valores estimados? Sé que hay muchos métodos para ajustar una curva suave, pero no estoy seguro de cuál sería el más apropiado para este tipo de curva y cómo lo escribirías R.

Me gusta loess()mucho suavizar:

x <- 1:10
y <- c(2,4,6,8,7,12,14,16,18,20)
lo <- loess(y~x)
plot(x,y)
lines(predict(lo), col='red', lwd=2)

El libro MASS de Venables y Ripley tiene una sección completa sobre suavizado que también cubre splines y polinomios, pero loess()es el favorito de todos.

Quizás smooth.spline es una opción. Puede establecer un parámetro de suavizado (normalmente entre 0 y 1) aquí

smoothingSpline = smooth.spline(x, y, spar=0.35)
plot(x,y)
lines(smoothingSpline)

también puede utilizar predecir en objetos smooth.spline. La función viene con la base R, consulte? Smooth.spline para obtener más detalles.

Para que sea REALMENTE suave ...

x <- 1:10
y <- c(2,4,6,8,7,8,14,16,18,20)
lo <- loess(y~x)
plot(x,y)
xl <- seq(min(x),max(x), (max(x) - min(x))/1000)
lines(xl, predict(lo,xl), col='red', lwd=2)

Este estilo interpola muchos puntos extra y le da una curva que es muy suave. También parece ser el enfoque que adopta ggplot. Si el nivel estándar de suavidad está bien, puede usarlo.

scatter.smooth(x, y)

la función qplot () en el paquete ggplot2 es muy simple de usar y proporciona una solución elegante que incluye bandas de confianza. Por ejemplo,

qplot(x,y, geom='smooth', span =0.5)

produce

LOESS es un muy buen enfoque, como dijo Dirk.

Otra opción es usar splines Bezier, que en algunos casos pueden funcionar mejor que LOESS si no tiene muchos puntos de datos.

Aquí encontrará un ejemplo: http://rosettacode.org/wiki/Cubic_bezier_curves#R

# x, y: the x and y coordinates of the hull points
# n: the number of points in the curve.
bezierCurve <- function(x, y, n=10)
    {
    outx <- NULL
    outy <- NULL

    i <- 1
    for (t in seq(0, 1, length.out=n))
        {
        b <- bez(x, y, t)
        outx[i] <- b$x
        outy[i] <- b$y

        i <- i+1
        }

    return (list(x=outx, y=outy))
    }

bez <- function(x, y, t)
    {
    outx <- 0
    outy <- 0
    n <- length(x)-1
    for (i in 0:n)
        {
        outx <- outx + choose(n, i)*((1-t)^(n-i))*t^i*x[i+1]
        outy <- outy + choose(n, i)*((1-t)^(n-i))*t^i*y[i+1]
        }

    return (list(x=outx, y=outy))
    }

# Example usage
x <- c(4,6,4,5,6,7)
y <- 1:6
plot(x, y, "o", pch=20)
points(bezierCurve(x,y,20), type="l", col="red")

Las otras respuestas son todos buenos enfoques. Sin embargo, hay algunas otras opciones en R que no se han mencionado, incluidas lowessy approx, que pueden brindar mejores ajustes o un rendimiento más rápido.

Las ventajas se demuestran más fácilmente con un conjunto de datos alternativo:

sigmoid <- function(x)
{
  y<-1/(1+exp(-.15*(x-100)))
  return(y)
}

dat<-data.frame(x=rnorm(5000)*30+100)
dat$y<-as.numeric(as.logical(round(sigmoid(dat$x)+rnorm(5000)*.3,0)))

Aquí están los datos superpuestos con la curva sigmoidea que los generó:

Este tipo de datos es común cuando se observa un comportamiento binario entre una población. Por ejemplo, esto podría ser un gráfico de si un cliente compró algo (un 1/0 binario en el eje y) versus la cantidad de tiempo que pasó en el sitio (eje x).

Se utilizan una gran cantidad de puntos para demostrar mejor las diferencias de rendimiento de estas funciones.

Smooth, spliney smooth.splinetodos producen un galimatías en un conjunto de datos como este con cualquier conjunto de parámetros que haya probado, tal vez debido a su tendencia a mapear en todos los puntos, lo que no funciona para datos ruidosos.

Las loess, lowessy approxlas funciones de todos los resultados utilizables producen, aunque apenas para approx. Este es el código para cada uno que usa parámetros ligeramente optimizados:

loessFit <- loess(y~x, dat, span = 0.6)
loessFit <- data.frame(x=loessFit$x,y=loessFit$fitted)
loessFit <- loessFit[order(loessFit$x),]

approxFit <- approx(dat,n = 15)

lowessFit <-data.frame(lowess(dat,f = .6,iter=1))

Y los resultados:

plot(dat,col='gray')
curve(sigmoid,0,200,add=TRUE,col='blue',)
lines(lowessFit,col='red')
lines(loessFit,col='green')
lines(approxFit,col='purple')
legend(150,.6,
       legend=c("Sigmoid","Loess","Lowess",'Approx'),
       lty=c(1,1),
       lwd=c(2.5,2.5),col=c("blue","green","red","purple"))

Como puede ver, lowessproduce un ajuste casi perfecto a la curva generadora original. Loessestá cerca, pero experimenta una extraña desviación en ambas colas.

Aunque su conjunto de datos será muy diferente, he descubierto que otros conjuntos de datos funcionan de manera similar, con ambos loessy lowesscapaces de producir buenos resultados. Las diferencias se vuelven más significativas cuando observa los puntos de referencia:

> microbenchmark::microbenchmark(loess(y~x, dat, span = 0.6),approx(dat,n = 20),lowess(dat,f = .6,iter=1),times=20)
Unit: milliseconds
                           expr        min         lq       mean     median        uq        max neval cld
  loess(y ~ x, dat, span = 0.6) 153.034810 154.450750 156.794257 156.004357 159.23183 163.117746    20   c
            approx(dat, n = 20)   1.297685   1.346773   1.689133   1.441823   1.86018   4.281735    20 a  
 lowess(dat, f = 0.6, iter = 1)   9.637583  10.085613  11.270911  11.350722  12.33046  12.495343    20  b 

Loesses extremadamente lento, tarda 100 veces más que approx. Lowessproduce mejores resultados que approx, sin dejar de correr bastante rápido (15 veces más rápido que loess).

Loess También se empantana cada vez más a medida que aumenta el número de puntos, volviéndose inutilizables alrededor de 50.000.

EDITAR: La investigación adicional muestra que loessofrece mejores ajustes para ciertos conjuntos de datos. Si se trata de un conjunto de datos pequeño o el rendimiento no es una consideración, pruebe ambas funciones y compare los resultados.

En ggplot2 puedes suavizar de varias formas, por ejemplo:

library(ggplot2)
ggplot(mtcars, aes(wt, mpg)) + geom_point() +
  geom_smooth(method = "gam", formula = y ~ poly(x, 2)) 
ggplot(mtcars, aes(wt, mpg)) + geom_point() +
  geom_smooth(method = "loess", span = 0.3, se = FALSE) 

No vi este método mostrado, así que si alguien más está buscando hacer esto, descubrí que la documentación de ggplot sugería una técnica para usar el gammétodo que producía resultados similares a los que se usaban loesscuando se trabajaba con conjuntos de datos pequeños.

library(ggplot2)
x <- 1:10
y <- c(2,4,6,8,7,8,14,16,18,20)

df <- data.frame(x,y)
r <- ggplot(df, aes(x = x, y = y)) + geom_smooth(method = "gam", formula = y ~ s(x, bs = "cs"))+geom_point()
r

Primero con el método loess y la fórmula automática Segundo con el método gam con la fórmula sugerida