Primero cargamos la base de datos desde la hoja de cálculo de Google

library("googlesheets4"); library(tidyverse); library(zoo); library(forecast); library(urca); library(tseries); library(tinytex)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.0      ✔ purrr   1.0.1 
## ✔ tibble  3.1.8      ✔ dplyr   1.0.10
## ✔ tidyr   1.3.0      ✔ stringr 1.5.0 
## ✔ readr   2.1.3      ✔ forcats 0.5.2 
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## 
## Attaching package: 'zoo'
## 
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## 
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
PORK <-read_sheet(
  'https://docs.google.com/spreadsheets/d/1hriCBeA-cjswm4_dxXWDZkxvffiikj3kxuSmsyu2efQ/edit#gid=1598108558',
  sheet = "s12",
  range = "s12!A6:B121",
  col_names = TRUE,
  col_types = NULL,
  na = "",
  trim_ws = TRUE,
  skip = 0,
  .name_repair = "unique"
)
## ! Using an auto-discovered, cached token.
##   To suppress this message, modify your code or options to clearly consent to
##   the use of a cached token.
##   See gargle's "Non-interactive auth" vignette for more details:
##   <]8;;https://gargle.r-lib.org/articles/non-interactive-auth.htmlhttps://gargle.r-lib.org/articles/non-interactive-auth.html]8;;>
## ℹ The googlesheets4 package is using a cached token for ']8;;mailto:majoa1204@gmail.commajoa1204@gmail.com]8;;'.
## ✔ Reading from "Series de tiempo - Econometria2".
## ✔ Range ''s12'!A6:B121'.

Ahora procedemos a modificar el formato de la fecha y a crear el logaritmo del precio promedio de cerdo en pie en Colombia ($ por kilo)

library("lubridate")
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
dataPORK <- PORK |>
  mutate(time = dmy(PORK$Fecha),
         precio = as.numeric(as.character(`Precio Promedio`))) |>
  select(time, precio) |>
  mutate(lprecio = log(precio))

Después de crear la nueva serie con mutate, usamos la funcion plot para observar el comportamiento de la variable

library(zoo)

{plot(dataPORK$time,dataPORK$lprecio, main = "",
     xlab = "Year quarter", ylab = "Precio de cerdo en pie",
     lwd = 0.4)
abline(lm(dataPORK$lprecio ~ dataPORK$time, data = dataPORK), col = "purple")}

Despues de modificar las variables, realizamos el grafico con el Log(precio) para analizar graficamente el comportamiento de la variable observada

ggplot(dataPORK) + geom_line(aes(time, log(precio), color = "Log(precio)"), linewidth = 0.7) + scale_color_manual(name = "", values = c("Log(precio)" = "darkviolet")) +labs(x="Fecha", y="Cerdo en pie") + theme(legend.text = element_text(size = 10), axis.title.x = element_text(size = 10), text = element_text(size=10))

Del gráfico es porsible observar que los datos son no estacionarios, con notables estacionalidad (cíclica) y tendencia no lineal

Análisis de Correlograma

ggAcf(dataPORK$lprecio, main="AC de LPrecio")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main`

ggPacf(dataPORK$lprecio, main="PAC de LPrecio")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main`

El correlograma nos indica que se trata de un proceso AR(2)

Posterior a los test y el análisis de correlogramas, realizamos los test de raíces unitarias y el test de Dickey-Fuller aumentado

 Para comenzar, se estima una regresión con una constante y una tendencia, adicionalmente, se selecciona AIC como criterio de información para seleccionar el número de rezagos óptimos

Test con de Dickey-Fuller con tendencia

DF1 <- ur.df(dataPORK$lprecio, type="trend", selectlags = "AIC")
summary(DF1)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression trend 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.155555 -0.016590  0.001651  0.017215  0.123433 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.7210698  0.2022801   3.565 0.000542 ***
## z.lag.1     -0.0864546  0.0242278  -3.568 0.000535 ***
## tt           0.0003903  0.0001373   2.843 0.005340 ** 
## z.diff.lag   0.6631310  0.0730085   9.083 5.27e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03394 on 109 degrees of freedom
## Multiple R-squared:  0.4453, Adjusted R-squared:  0.4301 
## F-statistic: 29.17 on 3 and 109 DF,  p-value: 6.35e-14
## 
## 
## Value of test-statistic is: -3.5684 4.5114 6.487 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2  6.22  4.75  4.07
## phi3  8.43  6.49  5.47
plot(DF1)

Para el test de Dickey-fuller con tendencia, los resultados arrojan que para un nivel de significancia de 90%, se rechaza la hipotesis nula de phi3 = (B1, B2, delta) = (B1, 0, o) - El rechazar phi3 puede implicar que delta es distinto de 0 y que sí hay una tendencia en el tiempo. En cuanto a phi2, este también se rechaza, pues el estadístico (=4.51) es mayor que el coeficiente obtenido al 10 (=4.07) - Rechazar phi2 implica que una, dos o todas estas tres variables son diferentes de cero

Test con de Dickey-Fuller con deriva

DF2 <- ur.df(dataPORK$lprecio, type="drift", selectlags = "AIC")
summary(DF2)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression drift 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.145390 -0.020733 -0.001816  0.014232  0.146937 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.33126    0.15343   2.159   0.0330 *  
## z.lag.1     -0.03836    0.01789  -2.144   0.0342 *  
## z.diff.lag   0.64822    0.07513   8.628 5.31e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03501 on 110 degrees of freedom
## Multiple R-squared:  0.4042, Adjusted R-squared:  0.3934 
## F-statistic: 37.31 on 2 and 110 DF,  p-value: 4.271e-13
## 
## 
## Value of test-statistic is: -2.1439 2.5613 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.46 -2.88 -2.57
## phi1  6.52  4.63  3.81
plot(DF2)

Este test estadístico de phi1 (H0: B1=delta=0) es de 2.56, menor a todos los valores críticos (6.52, 4.63, 3.81), por lo que no rechazo H0. Al no poder rechazar, se asume que tando B1 como delta son iguales a 0 y que 1) delta=0 por lo que hay una raíz unitaria presente y, 2) B1 = 0 por lo que no hay término de deriva.

Por otro lado, el tau2, muestra que -2.14 es menor que los valores críticos (-3.46, -2.88, -2.57), lo que indica que lprecio tiene una raíz unitaria

Test de Dickey-fuller sin tendencia y sin deriva

DF3 <- ur.df(dataPORK$lprecio, type="none", selectlags = "AIC")
summary(DF3)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.146996 -0.019355 -0.001069  0.017978  0.153729 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    0.0002625  0.0003930   0.668    0.506    
## z.diff.lag 0.6161187  0.0748466   8.232 3.95e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03558 on 111 degrees of freedom
## Multiple R-squared:  0.3886, Adjusted R-squared:  0.3776 
## F-statistic: 35.27 on 2 and 111 DF,  p-value: 1.386e-12
## 
## 
## Value of test-statistic is: 0.6679 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62
plot(DF3)

Para el test de Dickey-fuller sin deriva y sin tendencia, se tiene que la hipotesis nula es para delta = 0. El estadístico de prueba cae dentro de las 3 regiones por lo que no se puede rechazar la hipotesis nula - por esto mismo, se asume que los datos son de tipo caminata aleatoria, es decir, que hay una raíz unitaria presente

Probar serie diferenciada para ver si se logra estacionariedad

dataPORK <- dataPORK |> 
  mutate(dlprecio = c(NA,diff(dataPORK$lprecio)))

DF4 <- ur.df(dataPORK[!is.na(dataPORK$dlprecio),]$dlprecio, type="trend", selectlags = "AIC")
summary(DF4)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression trend 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.15418 -0.01518  0.00062  0.01357  0.10693 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.425e-03  6.141e-03  -0.232    0.817    
## z.lag.1     -5.662e-01  7.499e-02  -7.550 1.47e-11 ***
## tt           7.608e-05  9.346e-05   0.814    0.417    
## z.diff.lag   4.643e-01  8.524e-02   5.447 3.25e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03185 on 108 degrees of freedom
## Multiple R-squared:  0.3668, Adjusted R-squared:  0.3492 
## F-statistic: 20.85 on 3 and 108 DF,  p-value: 9.871e-11
## 
## 
## Value of test-statistic is: -7.5504 19.0095 28.5141 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2  6.22  4.75  4.07
## phi3  8.43  6.49  5.47
plot(DF4)

Los estadísticos caen por fuera de las 3 regiones para phi3, phi2 y tau3, por lo que se rechaza que es raíz unitaria

Realizar Test de Phillips-Perron y Test de Elliot-Rothenberg-Stock

# Test de Phillips-Perron
PP1 <- ur.pp(dataPORK$lprecio, type="Z-tau", model="trend", lags="short")
summary(PP1)
## 
## ################################## 
## # Phillips-Perron Unit Root Test # 
## ################################## 
## 
## Test regression with intercept and trend 
## 
## 
## Call:
## lm(formula = y ~ y.l1 + trend)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.181001 -0.023177  0.005149  0.023891  0.158625 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.3852458  0.2675690    1.44   0.1527    
## y.l1        0.9556918  0.0312038   30.63   <2e-16 ***
## trend       0.0003079  0.0001780    1.73   0.0864 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04459 on 111 degrees of freedom
## Multiple R-squared:  0.9473, Adjusted R-squared:  0.9464 
## F-statistic: 998.6 on 2 and 111 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic, type: Z-tau  is: -2.2156 
## 
##            aux. Z statistics
## Z-tau-mu             -0.1809
## Z-tau-beta            2.9048
## 
## Critical values for Z statistics: 
##                      1pct      5pct     10pct
## critical values -4.040722 -3.449402 -3.149694
plot(PP1)

# Test de Elliot-Rothenberg-Stock
ERS <- ur.ers(dataPORK$lprecio, type = "P-test", model = "trend", lag = 1)
summary(ERS)
## 
## ############################################### 
## # Elliot, Rothenberg and Stock Unit Root Test # 
## ############################################### 
## 
## Test of type P-test 
## detrending of series with intercept and trend 
## 
## Value of test-statistic is: 3.531 
## 
## Critical values of P-test are:
##                 1pct 5pct 10pct
## critical values 4.05 5.66  6.86

El test de Phillips-Perron plantea la hipótesis nula de que hay una raíz unitaria y la serie no es estacionaria, con la hipótesis alternativa de que no lo hay - Dado que el p-valor está por debajo de tamaño crítico, se rechaza la hipotesis nula

Finalmente, el estadístico del test de Elliot-Rothenberg-Stock es de 3.53, menor para todos los niveles de significancia (4.01 en el 1, 5.66 en el 5, 6.86 en el 10

MODELO ARIMA ESTACIONAL

library(urca); library(tidyverse); library(forecast); library(stats); library(lmtest)


ycerdo <- ts(na.omit(dataPORK$precio), start=2013, end=2022, frequency = 1)

ts.plot(ycerdo, main = "Log precio de cerdo en pie ($ por kilo)", xlab = "Año", col = "blue")
abline(h = mean(ycerdo), col = "orange", lwd = 2.5)

Primera diferencia Log precio de cerdo en pie en Colombia

ts.plot(diff(ycerdo), main = "Primera diferencia Log de cerdo en pie en Colombia", 
        xlab = "year") +
  abline(h = 0, col = "violet", lwd = 3 )

## integer(0)

Realizar correlograma de la serie para analizar el tipo de modelo

ggAcf(ycerdo, main = "Funcion AC log precio de cerdo en pie", ylab="AC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

ggPacf(ycerdo, main="Funcion PAC log precio de cerdo en pie", ylab="PAC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

A traves de los correlogramas podemos identificar un proceo AR(1), procederemos a realizar un proceso ARMA(1,0)

arma10 <- arima(ycerdo, order = c(1,0,0))
arma10
## 
## Call:
## arima(x = ycerdo, order = c(1, 0, 0))
## 
## Coefficients:
##          ar1  intercept
##       0.8917  4616.6873
## s.e.  0.1188   124.5477
## 
## sigma^2 estimated as 4612:  log likelihood = -57.16,  aic = 120.33
coeftest(arma10)
## 
## z test of coefficients:
## 
##             Estimate Std. Error z value  Pr(>|z|)    
## ar1          0.89174    0.11876  7.5087  5.97e-14 ***
## intercept 4616.68728  124.54770 37.0676 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

El arima no estacional (1,0,0) arroja que los parametros estimados son estadisticamente significativos y los valores satisfacen la condicion de estabilidad (|p| < 1)

Despues de esto, se revisan los residuales del modelo, el test Ljung-Box y el test Shapiro

#Residuales
res10 <- residuals(arma10)
ggAcf(res10, main="Funcion AC de los residuales", ylab="AC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

ggPacf(res10, main="Funcion PAC de los residuales", ylab="PAC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

#Realizar box test
Box.test(res10, lag = 1, type = "Ljung-Box")
## 
##  Box-Ljung test
## 
## data:  res10
## X-squared = 3.2495, df = 1, p-value = 0.07144
#Realizar Shapiro test
shapiro.test(res10)
## 
##  Shapiro-Wilk normality test
## 
## data:  res10
## W = 0.96769, p-value = 0.8686

El test de Ljung es utilizado básicamente para testear la ausencia de autocorrelacion serial hasta un rezago k especificado. El test determina si los errores son ruido blando o si hay algo más detrás de ellos - si las autocorrelaciones para los errores de los residuales son distintas a cero o no.

Lo anterior quiere decir que si la autocorrelacion de los residuales es muy pequeña, se puede decir que el modelo no exhibe una “falta significativa de ajuste” o que los residuales se distribuyen de manera independiente - Así pues, la H0 es que el modelo no muestra una falta de ajuste significativa, es decir, que el modelo está bien.

Este test arroja un estadistico de 3.24 y un p-valor=0.07, valor mayor que 0.05 por lo cual no se puede rechazar la hipótesis nula y concluir que los valores son independientes

Por otro lado, el test de Shapiro-Wilk prueba la normalidad de los residuales. Su H0 es que los residuales se distribuyen normalmente. En el presente caso se obtuvo un p-valor de 0.86 > 0.05, lo que indica que la distribución de la data no es significativamente diferente de una distribución normal, es decir, que distribuye normal y no se desvía de una distribución como esta

Se puede decir que los residuales no estan correlacionados y distribuyen normal

Para revisar los residuales más rápido de puede hacer uso de la funicon checkresiduals

checkresiduals(arma10)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,0) with non-zero mean
## Q* = 3.6727, df = 3, p-value = 0.299
## 
## Model df: 1.   Total lags used: 4

Se observa el comportamiento normal de los residuales

Se puede intentar realizar un modelo sobreparametrizado de ARMA(2,0)

arma200 <- arima(ycerdo, order=c(2,0,0))
arma200
## 
## Call:
## arima(x = ycerdo, order = c(2, 0, 0))
## 
## Coefficients:
##          ar1      ar2  intercept
##       1.4180  -0.6517  4597.6991
## s.e.  0.2252   0.2576    74.5258
## 
## sigma^2 estimated as 2831:  log likelihood = -55.15,  aic = 118.3
coeftest(arma200)
## 
## z test of coefficients:
## 
##             Estimate Std. Error z value  Pr(>|z|)    
## ar1          1.41796    0.22517  6.2974 3.027e-10 ***
## ar2         -0.65169    0.25761 -2.5297   0.01142 *  
## intercept 4597.69909   74.52583 61.6927 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Al comparar los resultados de los coeftest se observa que el coeficiente del primer rezago es similar al modelo ARMA(1,0), pero que el segundo rezago, a pesar de ser significativo, no lo es en la misma medida que el primero.

Para analizar el número de rezagos óptimos, se utilizan los criterior Akaike (AIC) y Bayesiano (BIC) - comparar arma10 y arma200

#AIC - AKAIKE
AIC(arma10);AIC(arma200)
## [1] 120.3291
## [1] 118.3045
#BIC - BAYESIANO
BIC(arma10); BIC(arma200)
## [1] 121.2369
## [1] 119.5148

Es claro que ambos criterios de informacion consideran al ARMA(2,0) como el mejor modelo - No obstante, para determinar con certeza cuál es el mejor modelo, se hace uso de un test likelihood-ratio - H0:modelo restringido (ARMA(1,0)) 

lrtestk <- as.numeric(2*(logLik(arma200)-logLik(arma10)))
pchisq(lrtestk, df=1, lower.tail = F)
## [1] 0.04484094

Con un p-valor de 0.04 se rechaza la hipotesis nula y se acepta el modelo ARMA (2,0) - Las mejoras SÍ son significativas para pasar de un modelo ARMA(1,0) a un modelo ARMA(2,0)

A continuacion se usara la función auto.arima para seleccionar el mejor modelo ARIMA para los datos estudiados

arma_pred <- auto.arima(ycerdo, stepwise = F, approximation = F)
arma_pred
## Series: ycerdo 
## ARIMA(0,1,0) 
## 
## sigma^2 = 4853:  log likelihood = -50.96
## AIC=103.92   AICc=104.49   BIC=104.12
forecast <- forecast::forecast(arma_pred, 10)
autoplot(forecast, main="Prediccion del precio de cerdo en pie en Colombia ($ por kilo)", col="red")

Con la funcion forecast, se es posible realizar una predicción del precio de cerdo en pie en Colombia ($ por kilo) para los próximos 10 años - Esta predicción se realiza con la serie de tiempo “ycerdo” y con el mejor modelo, seleccionado automáticamente por la funcion auto.arima el cual resulta ser un ARIMA(0,1,0)

Tal como lo arroja el auto.arima, el mejor modelo, con los menores criterios de informacion es un ARIMA(0,1,0)

Modelación ARIMA estacional extra

A continuación se describe la modelación de un ARIMA estacional con los datos de precio de cerdo en pie ($ por kilo). Se toman los datos para esta variable desde enero de 2013 a julio de 2022

#No se podria usar autoplot con dataPORK puesto que no esta funcion no soporta datos en data.frame - por esto se utliza "ycerdo" que anteriormente fue convertida en una seria de tiempo
ts.plot(ycerdo, main = "Log precio de cerdo en pie ($ por kilo)", xlab = "Año") +
  abline(h = mean(ycerdo), col = "red", lwd = 3)

## integer(0)

Se observa que los datos son no estacionarios, con notables estacionalidad (cíclica) y tendencia no lineal - Se toma una diferencia estacional en 12 meses para eliminar estacionalidad

library("tsibble"); library("dplyr")
## 
## Attaching package: 'tsibble'
## The following object is masked from 'package:lubridate':
## 
##     interval
## The following object is masked from 'package:zoo':
## 
##     index
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, union
dataPORK <- dataPORK |> mutate(precio.e12 = difference(dataPORK$precio, 12))
ts.plot(dataPORK$precio.e12)

Realizar correlogramas

forecast::ggAcf(dataPORK$precio.e12, main = "Funcion AC", ylab = "AC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

forecast::ggPacf(dataPORK$precio.e12, main = "Funcion PAC", ylab = "PAC")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main` and `ylab`

Los correlogramas con la diferencia estacional en 12 meses arrojan un modelo ARIMA(4,0) o ARMA (3,0). Adicionalmente, se observa que no hay estacionariedad - Se toman primeras diferencias y se vuelve a mirar el correlograma

dataPORK <- dataPORK |> mutate(dprecio.e12 = difference(dataPORK$precio.e12, 1))
ts.plot(dataPORK$dprecio.e12)

Dado que R no reconoce a la variable de precio como un vector parte de una serie de tiempo, incurrir en el autoarima es imposible pues lo reconoce como un dataframe. Es por esto que se mantiene la predicción de la funcion auto.arima para los proximos 10 años con respecto al precio de cerdo en pie ($ por kilo) en Colombia