Clase 8

Regresión polinomial dag = días de germinación pds = Profundidad de siembra

dag = sort(runif(35, min=6,max=12))
pds = rep(seq(0,12,2),each=5)

df = data.frame(dag=dag,pds=pds)
df

##          dag pds
## 1   6.227103   0
## 2   6.350429   0
## 3   6.387731   0
## 4   6.535562   0
## 5   6.802521   0
## 6   6.887542   2
## 7   7.056552   2
## 8   7.087829   2
## 9   7.375260   2
## 10  7.499431   2
## 11  7.689267   4
## 12  7.697036   4
## 13  7.718214   4
## 14  7.969930   4
## 15  8.197362   4
## 16  8.322907   6
## 17  8.375062   6
## 18  8.557371   6
## 19  8.918918   6
## 20  8.925418   6
## 21  9.089401   8
## 22  9.284820   8
## 23  9.403192   8
## 24 10.060801   8
## 25 10.282520   8
## 26 10.311275  10
## 27 10.434871  10
## 28 10.577286  10
## 29 10.871416  10
## 30 10.927412  10
## 31 11.419852  12
## 32 11.610662  12
## 33 11.619843  12
## 34 11.697315  12
## 35 11.843428  12

Modelo 1

mod1=lm(dag~pds, data=df)
summary(mod1)

## 
## Call:
## lm(formula = dag ~ pds, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.63207 -0.27478 -0.00806  0.24676  0.56105 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  6.26594    0.09568   65.49   <2e-16 ***
## pds          0.43194    0.01327   32.55   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.314 on 33 degrees of freedom
## Multiple R-squared:  0.9698, Adjusted R-squared:  0.9689 
## F-statistic:  1060 on 1 and 33 DF,  p-value: < 2.2e-16

Gráfico

plot(pds,dag,ylim = c(0,12))
abline(mod1)

Función

funpred=function(pds){
  if(0<=pds & pds <= 12){
    dag=6.91+0.43*pds
    return(dag)
  }else{"No está en el rango"}
}

funpred(2000)

## [1] "No está en el rango"

Modelo 2

mod2=lm(dag~poly(pds,2))
summary(mod2)

## 
## Call:
## lm(formula = dag ~ poly(pds, 2))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.50424 -0.18353 -0.05162  0.23502  0.68888 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    8.85759    0.04716  187.81  < 2e-16 ***
## poly(pds, 2)1 10.22159    0.27901   36.63  < 2e-16 ***
## poly(pds, 2)2  0.87323    0.27901    3.13  0.00372 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.279 on 32 degrees of freedom
## Multiple R-squared:  0.9769, Adjusted R-squared:  0.9754 
## F-statistic:   676 on 2 and 32 DF,  p-value: < 2.2e-16

Gráfica

plot(pds,dag,col='deepskyblue4',xlab='dps',main='Observed data')
lines(pds, fitted(mod2),col='firebrick1',lwd=3)
abline(mod1)

Función para la recta

funpred2 = function(pds){
  if(0<=pds & pds<=12){
    dag = 9.52+10.29*pds-1.56*(pds^2)
    return(dag)
  } else {"No está en el rango"}
}

Comparación entre modelos

library(Metrics)

## Warning: package 'Metrics' was built under R version 4.2.2

Metrics::rmse(dag, mod1$fitted.values)

## [1] 0.3048953

Metrics::rmse(dag, mod2$fitted.values) #mejor modelo

## [1] 0.2667855

#Akaike Information Criterion
AIC(mod1)

## [1] 22.18061

AIC(mod2)

## [1] 14.83398

AF = rnorm(n = 60,mean = 10,sd = 0.8)
fert = gl(3,20,60, c("D0", "D5", "D10"))

dm1 = data.frame(AF = AF, fert = fert)