Análisis y pronóstico de series de tiempo con R
Introducción
En este documento se explica lo relacionado al pronóstico de series de tiempo. Lo cual comprende:
Elaborando la serie de tiempo en R Análisis de la serie de tiempo Pronóstico de la serie de tiempo Estimación de errores de pronóstico Para el análisis de series de tiempo en R, existen una gran variedad de librerias. En este tutorial se emplea la libreria fpp2, la cual contiene contiene:
La libreria forecast, que contiene funciones de pronóstico La libreria ggplot2, que contiene funciones para gráficos Esta libreria se instala y se carga mediante install.packages(“fpp2”)y library(fpp2) respectivamente.
Elaborando la serie de tiempo en R Cargando librerias y datos Recordemos que se puede leer datos directamente desde una API, SQL u otros medios. Para este tutorial empleo la lectura de datos desde Excel. Entonces para ello emplearemos las siguientes funciones:
La función head() nos muestra los primeros datos de nuestra tabla de datos
library(fpp2)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo## ── Attaching packages ────────────────────────────────────────────── fpp2 2.5 ──## ✔ ggplot2 3.5.1 ✔ fma 2.5
## ✔ forecast 8.22.0 ✔ expsmooth 2.3## AirPassengers## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1949 112 118 132 129 121 135 148 148 136 119 104 118
## 1950 115 126 141 135 125 149 170 170 158 133 114 140
## 1951 145 150 178 163 172 178 199 199 184 162 146 166
## 1952 171 180 193 181 183 218 230 242 209 191 172 194
## 1953 196 196 236 235 229 243 264 272 237 211 180 201
## 1954 204 188 235 227 234 264 302 293 259 229 203 229
## 1955 242 233 267 269 270 315 364 347 312 274 237 278
## 1956 284 277 317 313 318 374 413 405 355 306 271 306
## 1957 315 301 356 348 355 422 465 467 404 347 305 336
## 1958 340 318 362 348 363 435 491 505 404 359 310 337
## 1959 360 342 406 396 420 472 548 559 463 407 362 405
## 1960 417 391 419 461 472 535 622 606 508 461 390 432Creando el objeto de serie de tiempo La serie de tiempo se debe almacenar en el objeto ts empleando la función ts(). El empleo de esta función es de la siguiente forma : serie <- ts(data, frequency= ,start=). Donde el primer agumento es un vector de datos. Luego los dos siguientes argumentos se indican de acuerdo a la siguiente tabla.
| Frecuencia | frequency= |
start= |
|---|---|---|
| Anual | 1 | 2000 |
| Trimestral | 4 | c(2000,2) |
| Mensual | 12 | c(2000,9) |
| Diario | 7 o 365.25 | 1 o c(2000,234) |
| Semanal | 52.18 | c(2000,23) |
| Horario | 24 o 168 o 8766 | 1 |
| Cada 30 min | 48 o 336 o 17532 | 1 |
En el siguiente fragmento de código definimos el objeto de la serie de tiempo y lo almacenaremos en data_serie.
file_path <- "C:/Users/Personal/OneDrive/Escritorio/SGA/Precipita.csv"
datavalencia<- read.csv(file_path)
datavalencia## DATOS
## 1 355.065
## 2 433.574
## 3 422.422
## 4 427.278
## 5 162.592
## 6 25.799
## 7 10.764
## 8 6.012
## 9 16.088
## 10 24.358
## 11 106.002
## 12 6.297
## 13 355.201
## 14 433.974
## 15 422.605
## 16 427.333
## 17 162.696
## 18 25.877
## 19 10.770
## 20 6.011
## 21 16.099
## 22 24.370
## 23 106.049
## 24 6.303
## 25 355.108
## 26 434.050
## 27 422.868
## 28 427.388
## 29 162.772
## 30 25.943
## 31 10.756
## 32 6.003
## 33 16.092
## 34 24.352
## 35 106.002
## 36 6.311
## 37 355.519
## 38 433.844
## 39 422.675
## 40 427.589
## 41 162.891
## 42 25.957
## 43 10.727
## 44 5.985
## 45 16.056
## 46 24.312
## 47 105.865
## 48 6.309
## 49 355.499
## 50 433.857
## 51 423.005
## 52 428.151
## 53 163.050
## 54 25.946
## 55 10.686
## 56 5.958
## 57 16.029
## 58 24.367
## 59 105.960
## 60 6.298
## 61 355.652
## 62 433.779
## 63 423.017
## 64 428.507
## 65 163.210
## 66 25.931
## 67 10.693
## 68 5.939
## 69 16.111
## 70 24.472
## 71 106.200
## 72 6.279
## 73 355.567
## 74 433.274
## 75 422.919
## 76 428.559
## 77 163.502
## 78 25.891
## 79 10.690
## 80 5.914
## 81 16.182
## 82 24.525
## 83 106.271
## 84 6.284
## 85 355.246
## 86 432.458
## 87 422.395
## 88 428.670
## 89 163.472
## 90 25.808
## 91 10.673
## 92 5.881
## 93 16.219
## 94 24.534
## 95 106.340
## 96 6.280
## 97 354.779
## 98 432.083
## 99 422.033
## 100 428.235
## 101 163.197
## 102 25.716
## 103 10.649
## 104 5.857
## 105 16.230
## 106 24.550
## 107 106.295
## 108 6.285
## 109 354.140
## 110 431.439
## 111 421.606
## 112 427.798
## 113 162.858
## 114 25.638
## 115 10.645
## 116 5.825
## 117 16.212
## 118 24.559
## 119 106.157
## 120 6.280
## 121 353.898
## 122 430.961
## 123 421.539
## 124 427.195
## 125 162.508
## 126 25.536
## 127 10.635
## 128 5.912
## 129 16.190
## 130 24.539
## 131 105.890
## 132 6.283
## 133 353.578
## 134 430.255
## 135 421.401
## 136 426.216
## 137 162.145
## 138 25.419
## 139 10.609
## 140 5.981
## 141 16.151
## 142 24.475
## 143 105.568
## 144 6.275
## 145 353.352
## 146 430.374
## 147 421.402
## 148 425.828
## 149 161.945
## 150 25.300
## 151 10.625
## 152 6.037
## 153 16.087
## 154 24.417
## 155 105.643
## 156 6.281
## 157 353.176
## 158 430.596
## 159 420.880
## 160 426.039
## 161 161.474
## 162 25.235
## 163 10.709
## 164 6.078
## 165 16.006
## 166 24.327
## 167 105.796
## 168 6.283
## 169 353.550
## 170 431.465
## 171 421.171
## 172 426.259
## 173 161.512
## 174 25.185
## 175 10.769
## 176 6.106
## 177 15.913
## 178 24.228
## 179 105.879
## 180 6.275
## 181 354.028
## 182 432.374
## 183 421.683
## 184 426.478
## 185 161.364
## 186 25.116
## 187 10.804
## 188 6.125
## 189 15.811
## 190 24.092
## 191 105.880
## 192 6.265
## 193 354.445
## 194 433.141
## 195 421.632
## 196 426.384
## 197 161.273
## 198 25.100
## 199 10.826
## 200 6.132
## 201 15.801
## 202 23.927
## 203 105.729
## 204 6.268
## 205 354.950
## 206 434.070
## 207 421.477
## 208 426.417
## 209 161.680
## 210 25.490
## 211 10.886
## 212 6.129
## 213 15.923
## 214 24.074
## 215 105.830
## 216 6.260
## 217 355.541
## 218 434.825
## 219 421.722
## 220 426.475
## 221 161.965
## 222 25.810
## 223 10.924
## 224 6.119
## 225 16.010
## 226 24.188
## 227 105.757
## 228 6.295
## 229 356.430
## 230 435.785
## 231 422.907
## 232 427.083
## 233 162.505
## 234 26.153
## 235 10.940
## 236 6.106
## 237 16.117
## 238 24.256
## 239 105.679
## 240 6.319
## 241 356.819
## 242 436.881
## 243 423.726
## 244 427.365
## 245 162.928
## 246 26.420
## 247 10.939
## 248 6.084
## 249 16.186
## 250 24.287
## 251 106.058
## 252 6.331
## 253 357.314
## 254 438.477
## 255 424.883
## 256 427.844
## 257 163.682
## 258 26.860
## 259 10.938
## 260 6.053
## 261 16.228
## 262 24.457
## 263 106.421
## 264 6.374
## 265 357.638
## 266 440.671
## 267 425.750
## 268 428.303
## 269 164.403
## 270 27.257
## 271 10.927
## 272 6.021
## 273 16.278
## 274 24.569
## 275 106.773
## 276 6.402
## 277 358.196
## 278 442.756
## 279 426.617
## 280 428.534
## 281 164.969
## 282 27.580
## 283 10.899
## 284 5.989
## 285 16.345
## 286 24.636
## 287 107.090
## 288 6.435
## 289 353.053
## 290 435.732
## 291 428.664
## 292 428.606
## 293 164.453
## 294 27.398
## 295 10.438
## 296 5.838
## 297 15.943
## 298 23.961
## 299 104.968
## 300 6.490
## 301 364.576
## 302 429.318
## 303 418.414
## 304 432.019
## 305 165.498
## 306 26.272
## 307 10.093
## 308 5.584
## 309 15.263
## 310 23.432
## 311 102.844
## 312 6.256
## 313 355.043
## 314 434.134
## 315 430.278
## 316 440.510
## 317 166.563
## 318 25.691
## 319 9.776
## 320 5.367
## 321 15.421
## 322 25.570
## 323 108.046
## 324 6.067
## 325 359.023
## 326 432.072
## 327 423.271
## 328 436.331
## 329 166.730
## 330 25.600
## 331 10.864
## 332 5.516
## 333 17.929
## 334 26.776
## 335 111.479
## 336 5.851
## 337 353.696
## 338 422.164
## 339 420.776
## 340 429.712
## 341 169.925
## 342 25.026
## 343 10.613
## 344 5.364
## 345 17.749
## 346 25.703
## 347 107.845
## 348 6.383
## 349 348.188
## 350 414.492
## 351 410.851
## 352 431.116
## 353 162.806
## 354 23.977
## 355 10.312
## 356 5.152
## 357 17.034
## 358 24.733
## 359 107.852
## 360 6.193
## 361 344.506
## 362 423.840
## 363 414.075
## 364 418.663
## 365 157.154
## 366 23.678
## 367 10.121
## 368 5.328
## 369 16.467
## 370 24.893
## 371 105.297
## 372 6.410
## 373 340.088
## 374 417.278
## 375 412.211
## 376 418.182
## 377 155.399
## 378 23.940
## 379 10.541
## 380 5.127
## 381 15.812
## 382 24.754
## 383 103.141
## 384 6.153
## 385 348.575
## 386 420.448
## 387 420.058
## 388 413.935
## 389 154.804
## 390 23.282
## 391 10.422
## 392 7.826
## 393 15.716
## 394 24.108
## 395 100.009
## 396 6.355
## 397 346.529
## 398 414.707
## 399 418.364
## 400 404.677
## 401 154.160
## 402 22.856
## 403 10.030
## 404 7.511
## 405 15.293
## 406 23.064
## 407 98.487
## 408 6.092
## 409 348.375
## 410 432.994
## 411 421.440
## 412 417.274
## 413 157.540
## 414 22.676
## 415 10.976
## 416 7.250
## 417 14.676
## 418 23.140
## 419 107.278
## 420 6.427
## 421 349.298
## 422 435.494
## 423 409.382
## 424 430.684
## 425 151.121
## 426 23.794
## 427 12.556
## 428 6.983
## 429 14.221
## 430 22.342
## 431 109.162
## 432 6.330
## 433 361.798
## 434 450.585
## 435 427.578
## 436 431.110
## 437 162.346
## 438 24.090
## 439 12.088
## 440 6.736
## 441 13.868
## 442 22.053
## 443 107.711
## 444 6.094
## 445 364.533
## 446 452.369
## 447 432.940
## 448 431.284
## 449 158.109
## 450 23.599
## 451 11.589
## 452 6.538
## 453 13.556
## 454 21.108
## 455 105.911
## 456 6.042
## 457 363.627
## 458 450.014
## 459 420.512
## 460 424.328
## 461 159.261
## 462 24.760
## 463 11.298
## 464 6.293
## 465 15.575
## 466 20.299
## 467 102.397
## 468 6.323
## 469 366.057
## 470 454.513
## 471 418.081
## 472 427.140
## 473 170.632
## 474 34.062
## 475 12.202
## 476 6.063
## 477 18.607
## 478 27.299
## 479 108.058
## 480 6.100
## 481 368.541
## 482 451.434
## 483 427.117
## 484 427.747
## 485 168.244
## 486 32.850
## 487 11.776
## 488 5.882
## 489 17.946
## 490 26.690
## 491 104.156
## 492 7.061
## 493 375.983
## 494 456.907
## 495 448.959
## 496 440.466
## 497 174.390
## 498 33.696
## 499 11.277
## 500 5.826
## 501 18.470
## 502 25.751
## 503 103.953
## 504 6.850
## 505 365.375
## 506 460.985
## 507 441.757
## 508 433.564
## 509 172.225
## 510 32.288
## 511 10.917
## 512 5.595
## 513 17.702
## 514 24.966
## 515 114.399
## 516 6.587
## 517 368.211
## 518 473.585
## 519 450.328
## 520 438.383
## 521 180.267
## 522 36.558
## 523 10.912
## 524 5.386
## 525 17.149
## 526 28.207
## 527 114.399
## 528 7.318
## 529 364.775
## 530 488.934
## 531 444.836
## 532 438.393
## 533 180.274
## 534 35.977
## 535 10.686
## 536 5.304
## 537 17.382
## 538 27.042
## 539 114.528
## 540 7.030
## 541 370.459
## 542 488.641
## 543 445.686
## 544 433.624
## 545 177.410
## 546 34.691
## 547 10.282
## 548 5.295
## 549 17.805
## 550 26.114
## 551 114.049
## 552 7.164
## 553 239.900
## 554 281.200
## 555 473.700
## 556 430.200
## 557 153.100
## 558 23.400
## 559 0.300
## 560 2.500
## 561 7.100
## 562 9.100
## 563 58.300
## 564 7.700
## 565 618.100
## 566 288.200
## 567 192.900
## 568 507.100
## 569 188.500
## 570 1.500
## 571 2.500
## 572 0.000
## 573 0.300
## 574 11.800
## 575 56.100
## 576 1.100
## 577 145.300
## 578 540.100
## 579 691.300
## 580 627.300
## 581 190.000
## 582 12.900
## 583 2.800
## 584 0.600
## 585 18.900
## 586 72.600
## 587 222.500
## 588 1.900
## 589 446.600
## 590 386.700
## 591 269.100
## 592 344.400
## 593 170.400
## 594 23.600
## 595 34.800
## 596 8.800
## 597 73.100
## 598 53.300
## 599 187.000
## 600 1.100
## 601 236.500
## 602 204.200
## 603 365.900
## 604 284.100
## 605 240.200
## 606 12.400
## 607 5.100
## 608 2.000
## 609 13.800
## 610 2.100
## 611 27.900
## 612 18.100
## 613 227.000
## 614 245.700
## 615 192.500
## 616 462.000
## 617 6.200
## 618 0.900
## 619 3.700
## 620 0.500
## 621 1.300
## 622 3.400
## 623 108.000
## 624 2.000
## 625 263.500
## 626 629.500
## 627 485.000
## 628 144.700
## 629 32.800
## 630 17.100
## 631 5.900
## 632 9.200
## 633 4.000
## 634 28.400
## 635 49.100
## 636 11.200
## 637 242.900
## 638 272.900
## 639 371.200
## 640 407.600
## 641 116.800
## 642 29.700
## 643 19.800
## 644 0.700
## 645 1.400
## 646 21.700
## 647 55.700
## 648 0.500
## 649 535.300
## 650 490.200
## 651 592.700
## 652 320.500
## 653 141.700
## 654 8.800
## 655 7.800
## 656 67.200
## 657 13.600
## 658 9.900
## 659 31.100
## 660 10.800
## 661 301.500
## 662 288.400
## 663 381.100
## 664 201.000
## 665 140.000
## 666 13.500
## 667 1.400
## 668 0.600
## 669 6.000
## 670 0.100
## 671 65.000
## 672 0.300
## 673 389.000
## 674 835.300
## 675 489.100
## 676 694.400
## 677 231.900
## 678 18.700
## 679 31.800
## 680 1.500
## 681 1.100
## 682 24.800
## 683 300.700
## 684 13.800
## 685 369.600
## 686 490.500
## 687 144.100
## 688 725.700
## 689 9.900
## 690 48.400
## 691 47.300
## 692 1.100
## 693 4.200
## 694 4.800
## 695 150.600
## 696 4.200
## 697 636.800
## 698 782.600
## 699 827.900
## 700 440.500
## 701 409.300
## 702 30.600
## 703 1.800
## 704 1.300
## 705 6.100
## 706 15.700
## 707 75.800
## 708 0.900
## 709 424.700
## 710 491.600
## 711 550.900
## 712 435.100
## 713 64.900
## 714 12.800
## 715 0.600
## 716 2.200
## 717 6.700
## 718 0.300
## 719 66.300
## 720 4.900
## 721 343.700
## 722 398.200
## 723 147.100
## 724 271.300
## 725 184.600
## 726 50.300
## 727 4.900
## 728 0.900
## 729 60.000
## 730 2.500
## 731 25.100
## 732 12.500
## 733 419.500
## 734 553.500
## 735 364.600
## 736 489.000
## 737 420.800
## 738 238.700
## 739 32.100
## 740 1.000
## 741 85.300
## 742 181.300
## 743 232.600
## 744 1.200
## 745 423.200
## 746 383.700
## 747 625.900
## 748 441.100
## 749 115.700
## 750 6.200
## 751 2.400
## 752 1.900
## 753 3.400
## 754 13.300
## 755 18.300
## 756 28.200
## 757 539.700
## 758 577.300
## 759 929.500
## 760 720.300
## 761 309.600
## 762 52.300
## 763 0.300
## 764 4.600
## 765 30.000
## 766 5.100
## 767 99.500
## 768 2.200
## 769 132.000
## 770 550.700
## 771 283.300
## 772 281.700
## 773 124.600
## 774 1.300
## 775 3.000
## 776 0.500
## 777 0.800
## 778 7.700
## 779 344.200
## 780 0.800
## 781 430.600
## 782 750.800
## 783 638.900
## 784 544.400
## 785 357.200
## 786 130.500
## 787 10.800
## 788 0.800
## 789 5.000
## 790 99.500
## 791 114.411
## 792 23.400
## 793 289.200
## 794 826.600
## 795 324.000
## 796 438.620
## 797 180.410
## 798 23.200
## 799 5.700
## 800 3.500
## 801 22.500
## 802 1.400
## 803 117.364
## 804 0.700
## 805 495.500
## 806 482.200
## 807 464.400
## 808 328.700
## 809 114.400
## 810 6.400
## 811 1.400
## 812 5.100
## 813 27.100
## 814 5.700
## 815 103.500
## 816 10.100View(datavalencia)
data_serie <- ts(datavalencia$DATOS, frequency=12, start=1946)
data_serie## Jan Feb Mar Apr May Jun Jul Aug Sep
## 1946 355.065 433.574 422.422 427.278 162.592 25.799 10.764 6.012 16.088
## 1947 355.201 433.974 422.605 427.333 162.696 25.877 10.770 6.011 16.099
## 1948 355.108 434.050 422.868 427.388 162.772 25.943 10.756 6.003 16.092
## 1949 355.519 433.844 422.675 427.589 162.891 25.957 10.727 5.985 16.056
## 1950 355.499 433.857 423.005 428.151 163.050 25.946 10.686 5.958 16.029
## 1951 355.652 433.779 423.017 428.507 163.210 25.931 10.693 5.939 16.111
## 1952 355.567 433.274 422.919 428.559 163.502 25.891 10.690 5.914 16.182
## 1953 355.246 432.458 422.395 428.670 163.472 25.808 10.673 5.881 16.219
## 1954 354.779 432.083 422.033 428.235 163.197 25.716 10.649 5.857 16.230
## 1955 354.140 431.439 421.606 427.798 162.858 25.638 10.645 5.825 16.212
## 1956 353.898 430.961 421.539 427.195 162.508 25.536 10.635 5.912 16.190
## 1957 353.578 430.255 421.401 426.216 162.145 25.419 10.609 5.981 16.151
## 1958 353.352 430.374 421.402 425.828 161.945 25.300 10.625 6.037 16.087
## 1959 353.176 430.596 420.880 426.039 161.474 25.235 10.709 6.078 16.006
## 1960 353.550 431.465 421.171 426.259 161.512 25.185 10.769 6.106 15.913
## 1961 354.028 432.374 421.683 426.478 161.364 25.116 10.804 6.125 15.811
## 1962 354.445 433.141 421.632 426.384 161.273 25.100 10.826 6.132 15.801
## 1963 354.950 434.070 421.477 426.417 161.680 25.490 10.886 6.129 15.923
## 1964 355.541 434.825 421.722 426.475 161.965 25.810 10.924 6.119 16.010
## 1965 356.430 435.785 422.907 427.083 162.505 26.153 10.940 6.106 16.117
## 1966 356.819 436.881 423.726 427.365 162.928 26.420 10.939 6.084 16.186
## 1967 357.314 438.477 424.883 427.844 163.682 26.860 10.938 6.053 16.228
## 1968 357.638 440.671 425.750 428.303 164.403 27.257 10.927 6.021 16.278
## 1969 358.196 442.756 426.617 428.534 164.969 27.580 10.899 5.989 16.345
## 1970 353.053 435.732 428.664 428.606 164.453 27.398 10.438 5.838 15.943
## 1971 364.576 429.318 418.414 432.019 165.498 26.272 10.093 5.584 15.263
## 1972 355.043 434.134 430.278 440.510 166.563 25.691 9.776 5.367 15.421
## 1973 359.023 432.072 423.271 436.331 166.730 25.600 10.864 5.516 17.929
## 1974 353.696 422.164 420.776 429.712 169.925 25.026 10.613 5.364 17.749
## 1975 348.188 414.492 410.851 431.116 162.806 23.977 10.312 5.152 17.034
## 1976 344.506 423.840 414.075 418.663 157.154 23.678 10.121 5.328 16.467
## 1977 340.088 417.278 412.211 418.182 155.399 23.940 10.541 5.127 15.812
## 1978 348.575 420.448 420.058 413.935 154.804 23.282 10.422 7.826 15.716
## 1979 346.529 414.707 418.364 404.677 154.160 22.856 10.030 7.511 15.293
## 1980 348.375 432.994 421.440 417.274 157.540 22.676 10.976 7.250 14.676
## 1981 349.298 435.494 409.382 430.684 151.121 23.794 12.556 6.983 14.221
## 1982 361.798 450.585 427.578 431.110 162.346 24.090 12.088 6.736 13.868
## 1983 364.533 452.369 432.940 431.284 158.109 23.599 11.589 6.538 13.556
## 1984 363.627 450.014 420.512 424.328 159.261 24.760 11.298 6.293 15.575
## 1985 366.057 454.513 418.081 427.140 170.632 34.062 12.202 6.063 18.607
## 1986 368.541 451.434 427.117 427.747 168.244 32.850 11.776 5.882 17.946
## 1987 375.983 456.907 448.959 440.466 174.390 33.696 11.277 5.826 18.470
## 1988 365.375 460.985 441.757 433.564 172.225 32.288 10.917 5.595 17.702
## 1989 368.211 473.585 450.328 438.383 180.267 36.558 10.912 5.386 17.149
## 1990 364.775 488.934 444.836 438.393 180.274 35.977 10.686 5.304 17.382
## 1991 370.459 488.641 445.686 433.624 177.410 34.691 10.282 5.295 17.805
## 1992 239.900 281.200 473.700 430.200 153.100 23.400 0.300 2.500 7.100
## 1993 618.100 288.200 192.900 507.100 188.500 1.500 2.500 0.000 0.300
## 1994 145.300 540.100 691.300 627.300 190.000 12.900 2.800 0.600 18.900
## 1995 446.600 386.700 269.100 344.400 170.400 23.600 34.800 8.800 73.100
## 1996 236.500 204.200 365.900 284.100 240.200 12.400 5.100 2.000 13.800
## 1997 227.000 245.700 192.500 462.000 6.200 0.900 3.700 0.500 1.300
## 1998 263.500 629.500 485.000 144.700 32.800 17.100 5.900 9.200 4.000
## 1999 242.900 272.900 371.200 407.600 116.800 29.700 19.800 0.700 1.400
## 2000 535.300 490.200 592.700 320.500 141.700 8.800 7.800 67.200 13.600
## 2001 301.500 288.400 381.100 201.000 140.000 13.500 1.400 0.600 6.000
## 2002 389.000 835.300 489.100 694.400 231.900 18.700 31.800 1.500 1.100
## 2003 369.600 490.500 144.100 725.700 9.900 48.400 47.300 1.100 4.200
## 2004 636.800 782.600 827.900 440.500 409.300 30.600 1.800 1.300 6.100
## 2005 424.700 491.600 550.900 435.100 64.900 12.800 0.600 2.200 6.700
## 2006 343.700 398.200 147.100 271.300 184.600 50.300 4.900 0.900 60.000
## 2007 419.500 553.500 364.600 489.000 420.800 238.700 32.100 1.000 85.300
## 2008 423.200 383.700 625.900 441.100 115.700 6.200 2.400 1.900 3.400
## 2009 539.700 577.300 929.500 720.300 309.600 52.300 0.300 4.600 30.000
## 2010 132.000 550.700 283.300 281.700 124.600 1.300 3.000 0.500 0.800
## 2011 430.600 750.800 638.900 544.400 357.200 130.500 10.800 0.800 5.000
## 2012 289.200 826.600 324.000 438.620 180.410 23.200 5.700 3.500 22.500
## 2013 495.500 482.200 464.400 328.700 114.400 6.400 1.400 5.100 27.100
## Oct Nov Dec
## 1946 24.358 106.002 6.297
## 1947 24.370 106.049 6.303
## 1948 24.352 106.002 6.311
## 1949 24.312 105.865 6.309
## 1950 24.367 105.960 6.298
## 1951 24.472 106.200 6.279
## 1952 24.525 106.271 6.284
## 1953 24.534 106.340 6.280
## 1954 24.550 106.295 6.285
## 1955 24.559 106.157 6.280
## 1956 24.539 105.890 6.283
## 1957 24.475 105.568 6.275
## 1958 24.417 105.643 6.281
## 1959 24.327 105.796 6.283
## 1960 24.228 105.879 6.275
## 1961 24.092 105.880 6.265
## 1962 23.927 105.729 6.268
## 1963 24.074 105.830 6.260
## 1964 24.188 105.757 6.295
## 1965 24.256 105.679 6.319
## 1966 24.287 106.058 6.331
## 1967 24.457 106.421 6.374
## 1968 24.569 106.773 6.402
## 1969 24.636 107.090 6.435
## 1970 23.961 104.968 6.490
## 1971 23.432 102.844 6.256
## 1972 25.570 108.046 6.067
## 1973 26.776 111.479 5.851
## 1974 25.703 107.845 6.383
## 1975 24.733 107.852 6.193
## 1976 24.893 105.297 6.410
## 1977 24.754 103.141 6.153
## 1978 24.108 100.009 6.355
## 1979 23.064 98.487 6.092
## 1980 23.140 107.278 6.427
## 1981 22.342 109.162 6.330
## 1982 22.053 107.711 6.094
## 1983 21.108 105.911 6.042
## 1984 20.299 102.397 6.323
## 1985 27.299 108.058 6.100
## 1986 26.690 104.156 7.061
## 1987 25.751 103.953 6.850
## 1988 24.966 114.399 6.587
## 1989 28.207 114.399 7.318
## 1990 27.042 114.528 7.030
## 1991 26.114 114.049 7.164
## 1992 9.100 58.300 7.700
## 1993 11.800 56.100 1.100
## 1994 72.600 222.500 1.900
## 1995 53.300 187.000 1.100
## 1996 2.100 27.900 18.100
## 1997 3.400 108.000 2.000
## 1998 28.400 49.100 11.200
## 1999 21.700 55.700 0.500
## 2000 9.900 31.100 10.800
## 2001 0.100 65.000 0.300
## 2002 24.800 300.700 13.800
## 2003 4.800 150.600 4.200
## 2004 15.700 75.800 0.900
## 2005 0.300 66.300 4.900
## 2006 2.500 25.100 12.500
## 2007 181.300 232.600 1.200
## 2008 13.300 18.300 28.200
## 2009 5.100 99.500 2.200
## 2010 7.700 344.200 0.800
## 2011 99.500 114.411 23.400
## 2012 1.400 117.364 0.700
## 2013 5.700 103.500 10.100plot(data_serie)
data_serie <- AirPassengers
data_serie## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1949 112 118 132 129 121 135 148 148 136 119 104 118
## 1950 115 126 141 135 125 149 170 170 158 133 114 140
## 1951 145 150 178 163 172 178 199 199 184 162 146 166
## 1952 171 180 193 181 183 218 230 242 209 191 172 194
## 1953 196 196 236 235 229 243 264 272 237 211 180 201
## 1954 204 188 235 227 234 264 302 293 259 229 203 229
## 1955 242 233 267 269 270 315 364 347 312 274 237 278
## 1956 284 277 317 313 318 374 413 405 355 306 271 306
## 1957 315 301 356 348 355 422 465 467 404 347 305 336
## 1958 340 318 362 348 363 435 491 505 404 359 310 337
## 1959 360 342 406 396 420 472 548 559 463 407 362 405
## 1960 417 391 419 461 472 535 622 606 508 461 390 432#Nota: Los resultados de esta guía corresponden a Data1 del archivo excel. Sugiero que tambien practiques con las otras series de tiempo como Data2 o Data3.
Análisis de la serie de tiempo Para el análisis de datos, esta libreria fpp2 nos permite emplear las siguientes funciones:
autoplot(), para graficar la serie de tiempo ggseasonplot(), para graficar la estacionalidad de una serie de tiempo ggsubseriesplot(), para graficar subseries ggAcf(), para graficar la autocorrelación decompose(), permite realizar una descomposición de estacionalidad y tendencia. Grafico de la serie de tiempo Para el gráfico simple de la serie de tiempo empleamos la funcion autoplot(). Donde solo se requiere ingresar el objeto de la serie de tiempo.
autoplot(data_serie)+
labs(title = "Serie de tiempo",
x = "Tiempo",
y = "Valor",
colour = "#00a0dc")+
theme_bw() 
Descomposición de la serie de tiempo Para la descomposición de la serie de tiempo se emplea la función decompose(). Donde se debe indicar el objeto de la serie de tiempo y el tipo de descomposición. Los tipos de descomposición que acepta esta función es additive y multiplicative.
Aditivo: Serie=T + S + I
Multiplicativo: Serie=T x S x I
Donde:
T: Tendencia
S: Estacionalidad
I: Irregular o error
# Descomposición de la serie de tiempo. Se almacena en el objeto fit
fit <- decompose(data_serie, type='additive')
#fit <- decompose(data_serie, type='multiplicative')
# Para graficar esta descomposición volvemos a emplear la funcion autoplot, pero con el objeto fit
autoplot(fit)+
labs(title = "Descomposicion de la serie de tiempo",
x = "Tiempo",
y = "Valor",
colour = "Gears")+
theme_bw()
library(highcharter)
hchart(stl(data_serie, s.window='periodic'))## Warning in create_yaxis(ntss, heights = heights, turnopposite = TRUE, title =
## list(text = NULL), : Deprecated function. Use the `create_axis` function.Grafico de la serie de tiempo con su tendencia El siguiente fragmento de código nos permite graficar la serie de tiempo con su tendencia. Notese que emplea el objeto fit en el cual guardamos previamente los valores de la descomposición. Nótese que se emplea la funcion trendcycle() para obtener los datos de tendencia del objeto fit.
autoplot(data_serie, series="Serie tiempo") +
autolayer(trendcycle(fit), series="Tendencia") +
labs(title = "Serie de tiempo",
x = "Tiempo",
y = "Valor"
) +
theme_bw()## Warning: Removed 12 rows containing missing values or values outside the scale range
## (`geom_line()`).
Grafico de estacionalidad Para realizar el gráfico de estacionalidad empleamos la función ggseasonplot. Donde el argumento es el objeto que contiene la serie de tiempo.
ggseasonplot(data_serie)
library(TSstudio)
ts_seasonal(data_serie, type = "all")ts_heatmap(data_serie)#Graficar la autocorrelacion
acf(data_serie, lag=35, col="2",main ="",
xlab="Lag",ylab="Valores")
#Graficar la autocorrelacion parcial
pacf(data_serie, lag=35, col="2",main ="",
xlab="Lag",ylab="Valores")
#Como hay confirmación de tendencia se hace la diferencia (1 lag)
plot(diff(data_serie),xlab="Tiempo",ylab="Valores")
acf(diff(data_serie), lag=35, col="2",main = "",
xlab="Lag", ylab="Valores")
pacf(diff(data_serie), lag=35,col="2",main = "",
xlab="Lag",ylab="Valores")
#Como la variancia se ve incrementa en el tiempo se debe sacar logaritmo natural de la serie
plot(diff(log(data_serie)),xlab="Tiempo",ylab="Valores")
acf(diff(log(data_serie)), lag=35, col="2",main = "",
xlab="Lag", ylab="Valores")
pacf(diff(log(data_serie)), lag=35,col="2",main = "",
xlab="Lag",ylab="Valores")
#Augmented Dickey-Fuller test para ver no estacionaridad
#Devuelve la probabilidad de ser no estacionaria
#Si es mayor que 0.05 la serie es no estacionaria (tiene tendencia)
#Si es menor que 0.05 es estacionaria (No hay tendencia)
# library(aTSA)
aTSA::adf.test(diff(log(data_serie)))## Augmented Dickey-Fuller Test
## alternative: stationary
##
## Type 1: no drift no trend
## lag ADF p.value
## [1,] 0 -9.61 0.01
## [2,] 1 -8.82 0.01
## [3,] 2 -7.63 0.01
## [4,] 3 -8.75 0.01
## [5,] 4 -6.79 0.01
## Type 2: with drift no trend
## lag ADF p.value
## [1,] 0 -9.63 0.01
## [2,] 1 -8.86 0.01
## [3,] 2 -7.71 0.01
## [4,] 3 -8.94 0.01
## [5,] 4 -6.98 0.01
## Type 3: with drift and trend
## lag ADF p.value
## [1,] 0 -9.60 0.01
## [2,] 1 -8.83 0.01
## [3,] 2 -7.69 0.01
## [4,] 3 -8.92 0.01
## [5,] 4 -6.95 0.01
## ----
## Note: in fact, p.value = 0.01 means p.value <= 0.01library(forecast)
#Ajuste automático de los coeficientes
ARIMAmodel <- auto.arima(log(data_serie))
ARIMAmodel## Series: log(data_serie)
## ARIMA(0,1,1)(0,1,1)[12]
##
## Coefficients:
## ma1 sma1
## -0.4018 -0.5569
## s.e. 0.0896 0.0731
##
## sigma^2 = 0.001371: log likelihood = 244.7
## AIC=-483.4 AICc=-483.21 BIC=-474.77#Ajuste automático de los coeficientes
ARIMAmodel <- auto.arima(diff(log(data_serie)))
ARIMAmodel## Series: diff(log(data_serie))
## ARIMA(0,0,1)(0,1,1)[12]
##
## Coefficients:
## ma1 sma1
## -0.4018 -0.5569
## s.e. 0.0896 0.0731
##
## sigma^2 = 0.001369: log likelihood = 244.7
## AIC=-483.39 AICc=-483.2 BIC=-474.77Pronóstico Métodos simples Para el pronóstico de series de tiempo mediante métodos básicos, la libreria fpp2 nos brinda las siguientes funciones:
naive(), metodo de naive simple ses(), exponential smoothing meanf(), media movil snaive(), metodo naive considerando estacionalidad El argumento a colocar en estas funciones es la serie de tiempo y el valor de h. Este valor de h es la cantidad de datos que deseamos pronosticar. Por ejemplo si deseamos pronosticar 12 datos, se debe indicar h=12.
Finalmente, para verificar el ajuste del método podemos emplear las siguientes funciones:
fitted(), obtiene un ajuste con la data historica checkresiduals(), permite analizar los residuales
# elaborando el método
m1 <- snaive(data_serie, h=96)
# graficando el pronóstico
autoplot(m1)
# verificando el ajuste del método
autoplot(m1)+autolayer(fitted(m1), series="Ajuste")## Warning: Removed 12 rows containing missing values or values outside the scale range
## (`geom_line()`).
# verificando los residuales
checkresiduals(m1)
##
## Ljung-Box test
##
## data: Residuals from Seasonal naive method
## Q* = 275.04, df = 24, p-value < 2.2e-16
##
## Model df: 0. Total lags used: 24Método regresión Para el pronóstico de series de tiempo mediante regresión, la libreria fpp2 nos brinda la función tslm(). Esta función la emplearemos para crear una regresión de la serie de tiempo con los datos de la descomposición estacional y/o tendencia. Entonces:
Si observamos solo tendencia usaremos : tslm(data_serie ~ trend) Si observamos solo estacionalidad usaremos : tslm(data_serie ~ season) Si observamos tendencia y estacionalidad usaremos : tslm(data_serie ~ trend + season) Luego con la función forecast realizamos el pronostico. El argumento a colocar en estas funcion es la regresión y el valor de h. Este valor de h es la cantidad de datos que deseamos pronosticar.
Finalmente, para verificar el ajuste del método podemos emplear las siguientes funciones:
fitted(), obtiene pronostico con la data historica checkresiduals(), permite analizar los residuales
library(tseries)
# elaborando la regresion
regresion <- tslm(data_serie ~ trend + season)
# elaborando el pronostico
m2 <- forecast(regresion)
# graficando el pronóstico
autoplot(m2)
# verificando el ajuste del método
autoplot(m2)+autolayer(fitted(m2), series="Ajuste")
# verificando los residuales
checkresiduals(m2)
##
## Ljung-Box test
##
## data: Residuals from Linear regression model
## Q* = 348.24, df = 24, p-value < 2.2e-16
##
## Model df: 0. Total lags used: 24Método holt winters Para el pronóstico de series de tiempo mediante holt winters, la libreria fpp2 nos brinda la función hw().
Los argumentos a colocar en esta funcion son:
La serie de tiempo El valor de h. Este valor de h es la cantidad de datos que deseamos pronosticar. El tipo de descomposición a usar para la estacionalidad. Los tipos de descomposición que acepta esta función es additive y multiplicative. Finalmente, para verificar el ajuste del método podemos emplear las siguientes funciones:
fitted(), obtiene pronostico con la data historica checkresiduals(), permite analizar los residuales
# elaborando el pronostico
m3 <- hw(data_serie, h=96, seasonal = 'multiplicative')
# graficando el pronóstico
autoplot(m3)
# verificando el ajuste del método
autoplot(m3)+autolayer(fitted(m3), series="Ajuste")
# verificando los residuales
checkresiduals(m3)
##
## Ljung-Box test
##
## data: Residuals from Holt-Winters' multiplicative method
## Q* = 40.151, df = 24, p-value = 0.0206
##
## Model df: 0. Total lags used: 24ARIMA Para el pronóstico de series de tiempo mediante ARIMA, la libreria fpp2 nos brinda la función auto.arima().
Primero crearemos un modelo ARIMA, para ello el argumento a colocar en esta funcion es la serie de tiempo. Considerar que esta función es solo una aproximación iterativa que busca los indices de AR y MA. Pues en determinadas series de tiempo podria no encontrar los indices adecuados para un modelo ARIMA. En ese caso lo adecuado es seguir la metodología de estimación de índices ARIMA. Esta metodología no esta cubierta en esta guia.
Luego con la función forecast realizamos el pronostico. El argumento a colocar en estas funcion es el modelo ARIMA y el valor de h. Este valor de h es la cantidad de datos que deseamos pronosticar.
Finalmente, para verificar el ajuste del método podemos emplear las siguientes funciones:
fitted(), obtiene pronostico con la data historica checkresiduals(), permite analizar los residuales
# elaborando el modelo ARIMA
modelo_arima <- auto.arima(data_serie)
# elaborando el pronostico
m4 <- forecast(modelo_arima, h=96)
# graficando el pronóstico
autoplot(m4)
# verificando el ajuste del método
autoplot(m4)+autolayer(fitted(m4), series="Ajuste")
# verificando los residuales
checkresiduals(m4)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,1)(0,1,0)[12]
## Q* = 37.784, df = 21, p-value = 0.01366
##
## Model df: 3. Total lags used: 24Red neuronal Para el pronóstico de series de tiempo mediante una red neuronal, la libreria fpp2 nos brinda la función nnetar().
Primero crearemos un modelo de red neuronal (neural network), para ello el argumento a colocar en esta funcion es la serie de tiempo.
Luego con la función forecast realizamos el pronostico. El argumento a colocar en estas funcion es el modelo de red neuronal y el valor de h. Este valor de h es la cantidad de datos que deseamos pronosticar.
Finalmente, para verificar el ajuste del método podemos emplear las siguientes funciones:
fitted(), obtiene pronostico con la data historica checkresiduals(), permite analizar los residuales
# elaborando el modelo de red neuronal
neural_network <- nnetar(data_serie)
# elaborando el pronostico
m5 <- forecast(neural_network, h=96)
# graficando el pronóstico
autoplot(m5)
# verificando el ajuste del método
autoplot(m5)+autolayer(fitted(m5), series="Ajuste")## Warning: Removed 12 rows containing missing values or values outside the scale range
## (`geom_line()`).
# verificando los residuales
checkresiduals(m5)
##
## Ljung-Box test
##
## data: Residuals from NNAR(1,1,2)[12]
## Q* = 193.67, df = 24, p-value < 2.2e-16
##
## Model df: 0. Total lags used: 24Estimación de error Para estimar los errores de pronóstico, se debe realizar con los valores ocurridos o reales. Es decir este error se mide tiempo despues de haber realizado el pronóstico.
Supongamos que los valores reales ocurridos para la Data1 son los siguientes:
Mes Valor Real
Enero 19 13487 Mayo Febrero 19 12776 Junio Marzo 19 13812 Julio Abril 19
13032 Agosto Entonces, para poder comparar nuestros datos sera necesario
almacenarlo en un objeto de serie de tiempo. Entonces:
real <- c(622, 606, 508, 461, 390, 432,462, 471)
data_real <- ts(real, frequency=12,start=1961)La función accuracy() determina los errores de pronostico, para ello es necesario:
Indicar el modelo de pronostico, donde estará el pronóstico Indicar datos reales, donde estará los valores reales Entonces lo que realizará esta función es comparar el pronóstico de los siguientes 8 datos y el valor real. Pues estamos considerando que ya pasaron 8 meses, y nos encontramos en la etapa de evaluar el error de pronóstico de nuestros modelos.
# modelo en base a métodos simples
accuracy(m1,data_real)## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 31.77273 36.31574 32.0303 11.124393 11.24871 1.000000 0.7464603
## Test set 3.62500 140.23596 123.6250 -2.775782 24.26489 3.859626 0.6988219
## Theil's U
## Training set NA
## Test set 2.475277# modelo en base a regresion lineal
accuracy(m2,data_real)## ME RMSE MAE MPE MAPE MASE
## Training set 1.999781e-16 25.11363 19.77364 0.6039558 8.591515 0.6173418
## Test set -9.952652e-01 106.24590 92.63684 -3.1086876 18.289105 2.8921624
## ACF1 Theil's U
## Training set 0.7631551 NA
## Test set 0.6491973 1.800463# modelo en base a holt winters
accuracy(m3,data_real)## ME RMSE MAE MPE MAPE MASE
## Training set 1.256973 10.63256 7.790649 0.2182707 2.914411 0.2432275
## Test set -35.415680 150.14883 136.598759 -10.9560051 27.796830 4.2646727
## ACF1 Theil's U
## Training set 0.2135914 NA
## Test set 0.6873988 3.013635# modelo en base a ARIMA
accuracy(m4,data_real)## ME RMSE MAE MPE MAPE MASE
## Training set 1.34230 10.84619 7.86754 0.420698 2.800458 0.245628
## Test set -27.27284 144.07536 132.46579 -9.203000 26.842386 4.135640
## ACF1 Theil's U
## Training set -0.00124847 NA
## Test set 0.69852717 2.839455# modelo en base a red neuronal
accuracy(m5,data_real)## ME RMSE MAE MPE MAPE MASE
## Training set -0.003222003 14.4341 11.56672 -0.298241 4.268011 0.3611179
## Test set -47.394919155 157.0424 143.84338 -13.475238 29.603631 4.4908528
## ACF1 Theil's U
## Training set 0.5539169 NA
## Test set 0.7017837 3.272749Entonces a partir de estos errores de pronóstico podemos determinar cual es el modelo mas adecuado. En la práctica se suelen evaluar varios modelos e incluso tomar como pronóstico un valor medio del resultados de dos o más modelos.