Clase 9

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.195180   0
## 2   6.461419   0
## 3   6.484119   0
## 4   6.503660   0
## 5   6.524377   0
## 6   6.599206   2
## 7   6.610643   2
## 8   6.705634   2
## 9   6.732643   2
## 10  6.867840   2
## 11  6.949977   4
## 12  7.003800   4
## 13  7.086481   4
## 14  7.115751   4
## 15  7.250028   4
## 16  7.314165   6
## 17  7.316944   6
## 18  7.375065   6
## 19  8.101912   6
## 20  8.227447   6
## 21  8.325886   8
## 22  8.392321   8
## 23  9.210441   8
## 24  9.306890   8
## 25  9.660385   8
## 26  9.796211  10
## 27 10.360651  10
## 28 10.490951  10
## 29 10.640237  10
## 30 10.904729  10
## 31 11.033765  12
## 32 11.422569  12
## 33 11.589275  12
## 34 11.923329  12
## 35 11.979894  12

Modelo 1

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

## 
## Call:
## lm(formula = dag ~ pds, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.09909 -0.40078  0.00661  0.44797  0.90556 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  5.75217    0.18346   31.35   <2e-16 ***
## pds          0.44351    0.02544   17.43   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6021 on 33 degrees of freedom
## Multiple R-squared:  0.902,  Adjusted R-squared:  0.8991 
## F-statistic: 303.9 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.71128 -0.20076  0.06123  0.22694  0.76253 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    8.41325    0.06271  134.16  < 2e-16 ***
## poly(pds, 2)1 10.49544    0.37099   28.29  < 2e-16 ***
## poly(pds, 2)2  2.74911    0.37099    7.41 1.98e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.371 on 32 degrees of freedom
## Multiple R-squared:  0.9639, Adjusted R-squared:  0.9617 
## F-statistic: 427.6 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.5846075

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

## [1] 0.354732

#Akaike Information Criterion
AIC(mod1)

## [1] 67.74868

AIC(mod2)

## [1] 34.7782

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)