##EtudedelaconcentrationdeCO2dansl'atmosphère depuis 1958
Fichier télécharger le 18/02/2021
Atmospheric CO2 concentrations (ppm) derived from in situ air measurements at Mauna Loa, Observatory, Hawaii: Latitude 19.5°N Longitude 155.6°W Elevation 3397m.
Source: R. F. Keeling, S. J. Walker, S. C. Piper and A. F. Bollenbacher. Scripps CO2 Program ( http://scrippsco2.ucsd.edu ). Scripps Institution of Oceanography (SIO). University of California. University of California. La Jolla, California USA 92093-0244.
Données hebdomadaires.
Les analyses effectuées pour décomposer la série temporelle ont été réalisées à partir du livre *Numerical Ecology* écrit par Pierre et Louis Lengendre et édité par Elsevier en 1998.
`url_data_co2_web` = lien pour télécharger les données.
### 1. Réalisez un graphique qui vous montrera une oscillation périodique superposée à une évolution systématique plus lente.
Après un examen visiuel du fichier csv, l'onconstatequelesdonnéescommencentàlaligne45.L'on peut donc "skipper" les 44 premières lignes. La première colonne correspond aux dates et la seconde aux concentration de CO2 en ppm.
On vérifie le type de données dans chaque colonne.
```{r}
class(data_co2$dates)
class(data_co2$CO2)
```
Il faut convertir le vecteur *dates* au format date avec la fonction `as.Date`. Les dates renseignées étant déjà au bon format, cela ne devrait pas poser de problème. Les résultats sont placés dans un vecteur *dates2*:
```{r}
data_co2$dates2 = as.Date(data_co2$dates, format = "%Y-%m-%d")
```
On vérifie que cela a fonctionné:
```{r}
class(data_co2$dates2)
```
Y a-t-il des points manquants dans nos données ?
```{r}
na_records = apply(data_co2, 1, function (x) any(is.na(x)))
data_co2[na_records,]
```
On peut à présent représenter graphiquement la relation entre la concentration en CO2 et le temps:
```{r}
plot(data_co2$CO2 ~ data_co2$dates2, type = "l",
col = "blue", xlab = "Time", ylab = "CO2 in ppm")
```
On zoom sur la fin des données pour mieux voir la périodicité.
```{r}
plot(data_co2$CO2[c(3000:3208)] ~ data_co2$dates2[c(3000:3208)], type = "l",
col = "blue", xlab = "Time", ylab = "CO2 in ppm")
```
On peut donc constater la présence d'unepériodicité,quisembleàpremièrevueêtreannuelle,etd'une tendance générale à l'augmentation.
###2.Séparezcesdeuxphénomènes.Caractérisezl'oscillation périodique. Proposez un modèle simple de la contribution lente, estimez ses paramètres et tentez une extrapolation jusqu'à2025(danslebutdepouvoirvaliderlemodèlepardesobservationsfutures).
Pourdécomposercettesérie,deuxmodèlesdedécompositions'offrent à nous: un modèle additif ou multiplicatif. La périodicité étant annuelle, nous allons estimer les moyennes et les écarts types annuels. Si aucune relation n'apparaitentrecesdeuxparamètres,onappliqueraunmodèleadditifpourdécomposerlasérie,sinononappliqueraunmodèlemultiplicatif.
Onchargelepackage`lubridate`quivanouspermettred'utiliser la fonction `year` que nous allons utiliser pour isoler les concentrations en CO2 qui appartiennent toutes à la même année.
```{r include=FALSE}
library(lubridate)
```
On créé deux vecteurs vides, `moyenne_ann` et `et_ann` pour stocker la concentration en CO2 moyenne de chaque année et l'écarttypeassocié
Oneffectueensuiteunebouclede1958(annéeminimale(`min(year(data_co2$dates2))`))à2021(annéemaximale(`max(year(data_co2$dates2))`)).Ondéfinitunconteur`a`égalà1dontlavaleurchangeraàchaqueitérationpournouspermettredestockerlesmoyennesetlesécartstypesàlabonnepositiondansnosvecteursvides`moyenne_ann`et`et_ann`.Ondéfinitégalementunvecteur*logique*`test`quivanouspermettred'isoler dans la matrice `data_co2$CO2` les concentration en CO2 correspondant à l'annéedel'itérations.
```{r}
a=1
for (i in min(year(data_co2$dates2)):max(year(data_co2$dates2))){
Et calculer le coefficient de corrélation de Pearson `r` entre les deux paramètres avec la fonction `cor`.
```{r}
r = cor(moyenne_ann, et_ann, method = "pearson")
r
```
Il n'yapasderelationentrelesdeuxvariables,onpeutdoncutiliserunmodèleadditifpourdécomposerlasérie.
Pourextrairelatendancedesdonnées,onappliqueunemoyennemobilesimpled'ordre n = 52 (nombre de semaines dans une année) avec la fonction `movavg` du package `pracma`
```{r}
library(pracma)
trend = movavg(data_co2$CO2, 52, type = "s")
```
On plot le résultat
```{r}
plot(trend ~ data_co2$dates2, type = "l", col = "red", xlab = "", ylab = "", axes = FALSE)
par(new = T)
plot(data_co2$CO2 ~ data_co2$dates2, type = "l", col = "blue",
On extrait ensuite la tendance des données avec la formule suivante $Y_{i Detrend} = Y_{i} - \hat{Y_{i}}$, avec $Y_{i}$ les concentrations en CO2 et $\hat{Y_{i}}$ la tendance estimé avec la moyenne mobile.
```{r}
data_co2$CO2dt = data_co2$CO2 - trend
```
On plot le résultat pour voir si cela a marché.
```{r}
plot(data_co2$CO2dt ~ data_co2$dates2, type = "l", col = "blue",
xlab = "Time", ylab = "Detrend CO2 in ppm")
```
Pour caractériser l'oscillationpériodique,nousallonseffectuéunemesured'autocorrelation avec la fonction `acf` sur les données détentancées:
```{r}
auto_cor = acf(data_co2$CO2dt, type = "correlation", lag.max = 100)
```
On observe un décalage d'environs52-53pasdetemps,cequicorrespondàladuréed'une année en semaine. Nous avons donc une périodicité annuelle.
On va ensuite ajuster un polynome de degré 3 avec la fonction `polyfit` du package `pracma`.
```{r}
model = polyfit(seq(1,3208), data_co2$CO2, 3)
```
Voici les paramètres du modèle:
```{r}
model
```
Et sa représentation graphique.
```{r}
output = polyval(model, seq(1,3208))
plot(output ~ data_co2$dates2, type = "l", col = "red", axes = FALSE,
xlab = "", ylab = "", lwd = 2)
par(new = T)
plot(data_co2$CO2 ~ data_co2$dates2, type = "l", col = "blue",
xlab = "Time", ylab = "CO2 in ppm")
legend(1, 400, legend = c("model", "raw data"),
col = c("red", "blue"), lty = 1, cex = 1)
```
Malheureusement, nous ne sommes pas parvenu à ajuster le polynome en considérant les dates en tant que tel.
La dernière valeur de CO2 enregistrée correspond à la date:
