Ajustar Curvas a los datos en R

Paso 1: Cargamos los datos que queremos fitear.

En este caso tenemos la variable Counts que depende de la variable Time.

dataToFit <- tibble(Time = c( 1, 2, 4, 6, 8, 9, 10, 11, 12, 13, 
14, 15, 16, 18, 19, 20, 21, 22, 24, 25, 26, 27, 28, 29, 30,31), 
Counts = c(126.6, 101.8, 71.6, 101.6, 68.1, 62.9, 45.5, 41.9, 
46.3, 34.1, 38.2, 41.7, 24.7, 41.5, 36.6, 19.6, 
22.8, 29.6, 23.5, 15.3, 13.4, 26.8, 9.8, 18.8, 25.9, 19.3))
ggplot(dataToFit, aes(x=Time,y=Counts))+
  geom_point()

Podemos ver que estos datos siguen una tendencia exponencial. Por lo tanto, trataremos de ajustar una función de estas característica a los datos.

Paso 2: Creamos la función y calculamos los parámetros que mejor se ajustan.

Crearemos una función exponencial con la siguiente forma \(f(x)=b.x^{c}\)

exponencial <- function(x,b,c){return(b*x^(c))}

Podemos ver que en esta función tenemos dos constantes que desconocemos: b y c. Para encontrar esas constantes que mejor describen la relación entre las variables dependiente e independiente de nuestro problema. Para eso usaremos la función nls que recibe los siguientes parámetros:

• data: dataframe que contiene los datos
• formula: mapeo de las variables a relacionar con la siguiente forma y~f(x,params) donde
  • y: variable dependiente que es una columna del dataset
  • f: función definida previamente que se quiere ajustar a los datos
  • x: variable independiente que es una columna del dataset
  • params: conjunto de parámetros que se quieren ajustar a los datos
• start: lista de valores iniciales para los parámetros que se quieren encontrar

En nuestro caso particular los valores de estos parámetros serían

• data: dataToFit
• formula:
  • y: Counts
  • f: exponencial (definida en el chunk anterior)
  • x: Time
  • params: b y c
• start: b=6 y c=0.5

ajuste <- nls(data=dataToFit, formula = Counts~exponencial(Time,b,c),start=list(b=6,c=0.5))

Si el parámetro start está mal definido podría causar que la función nls nunca encuentre parámetros correctos. Por lo tanto un buen truco para definirlos es agarrar un datapoint y encontrar los valores para los parámetros. En nuestro caso si \(Time = 16\), \(Counts \approx 24\) . Por lo tanto usando la función exponencial podríamos decir que \(24=6.16^{0.5}\) y de esa manera concluir que b=6 y c=0.5 podrían ser buenos valores iniciales.

Paso 3: Recuperar parámetros y extrapolar valores

La variable ajuste contiene los resultados de la función nls. Podemos acceder a los parámetros encontrados de la siguiente manera.

(b_ajustado <- ajuste$m$getPars()["b"])

##        b 
## 139.8903

(c_ajustado <- ajuste$m$getPars()["c"])

##          c 
## -0.4900971

Ahora podemos usar esos valores ajustados para interpolar y extrapolar valores de la función que no conocíamos.

(fitData <- tibble(Time=1:50,fit = exponencial(Time,b_ajustado,c_ajustado)))

En el siguiente gráfico podemos ver qué tan bien se ajusta la curva a los datos originales.

ggplot()+
  geom_point(dataToFit,mapping=aes(x=Time,y=Counts)) +
  geom_line(fitData,mapping=aes(x=Time,y=fit))