Background
In a paper from Epstein (1991): On obtaining daily climatological values from monthly means, the formulation and an algorithm for calculating Fourier interpolation for periodical and even-spaced values are given.
In the paper, the goal is to obtain daily values from monthly means by interpolation.
In short, it is assumed that unknown daily values can be represented by the sum of harmonic components: $$ y(t) = a_{0} + \sum_{j}\left[a_{j}\,\cos(2\pi jt/12)+b_{j}\,\sin(2\pi jt/12)\right] $$ In the paper $t$ (time) is expressed in months.
After some derviation, it is shown that the terms can be calculated by: $$ \begin{align} a_{0} &= \sum_{T}Y_{T}/12 \\ a_{j} &= \left[ (\pi j/12)/\sin(\pi j/12)\right] \times \sum_{T}\left[Y_{T}\,\cos(2\pi jT/12)/6 \right]~~~~~~~j=1,\ldots, 5 \\ b_{j} &= \left[ (\pi j/12)/\sin(\pi j/12)\right] \times \sum_{T}\left[Y_{T}\,\sin(2\pi jT/12)/6 \right]~~~~~~~j=1,\ldots, 5 \\ a_{6} &= \left[ (\pi j/12)/\sin(\pi j/12)\right]\times \sum_{T}\left[Y_{T}\cos(\pi T)/12\right] \\ b_{6} &= 0 \end{align} $$ Where $Y_{T}$ denote the monthly means and $T$ the month.
Harzallah (1995) summarizes this aproach as follows: "The interpolation is carried out by adding zeros to the spectral coefficients of data and by performing an inverse Fourier transform to the resulting extended coefficients. The method is equivalent to applying a rectangular filter to Fourier coefficients."
Questions
My goal is to use the above methodology for interpolation of weekly means to obtain daily data (see my previous question). In summary, I have 835 weekly means of count data (see the example dataset at the bottom of the question). There are quite a few things that I don't understand before I can apply the approach outlined above:
- How would the formulas have to be changed for my situation (weekly instead of monthly values)?
- How could the time $t$ be expressed? I assumed $t/835$ (or $t/n$ with $n$ data points in general), is that correct?
- Why does the author calculate 7 terms (i.e. $0\leq j \leq 6$)? How many terms would I have to consider?
- I understand that the question can probably be solved by using a regression approach and using the predictions for interpolation (thanks to Nick). Still, some things are unclear to me: How many terms of harmonics should be included in the regression? And what period should I take? How can the regression be done to ensure that the weekly means are preserved (as I don't want an exact harmonic fit to the data)?
Using the regression approach (which is also explained in this paper), I managed to get an exact harmonic fit to the data (the $j$ in my example would run through $1, \ldots, 417$, so I fitted 417 terms). How can this approach be modified -$~$if possible$~$- to achieve the conservation of the weekly means? Maybe by applying correction factors to each regression term?
The plot of the exact harmonic fit is:
EDIT
Using the signal package and the interp1
function, here's what I've managed to do using the example data set from below (many thanks to @noumenal). I use q=7
as we have weekly data:
# Set up the time scale
daily.ts <- seq(from=as.Date("1995-01-01"), to=as.Date("2010-12-31"), by="day")
# Set up data frame
ts.frame <- data.frame(daily.ts=daily.ts, wdayno=as.POSIXlt(daily.ts)$wday,
yearday = 1:5844,
no.influ.cases=NA)
# Add the data from the example dataset called "my.dat"
ts.frame$no.influ.cases[ts.frame$wdayno==3] <- my.dat$case
# Interpolation
case.interp1 <- interp1(x=ts.frame$yearday[!is.na(ts.frame$no.influ.case)],y=(ts.frame$no.influ.cases[!is.na(ts.frame$no.influ.case)]),xi=ts.frame$yearday, method = c("cubic"))
# Plot subset for better interpretation
par(bg="white", cex=1.2, las=1)
plot((ts.frame$no.influ.cases)~ts.frame$yearday, pch=20,
col=grey(0.4),
cex=1, las=1,xlim=c(0,400), xlab="Day", ylab="Influenza cases")
lines(case.interp1, col="steelblue", lwd=1)
There are two issues here:
- The curve seem to fit "too good": it goes through every point
- The weekly means are not conserved
Example dataset
structure(list(date = structure(c(9134, 9141, 9148, 9155, 9162,
9169, 9176, 9183, 9190, 9197, 9204, 9211, 9218, 9225, 9232, 9239,
9246, 9253, 9260, 9267, 9274, 9281, 9288, 9295, 9302, 9309, 9316,
9323, 9330, 9337, 9344, 9351, 9358, 9365, 9372, 9379, 9386, 9393,
9400, 9407, 9414, 9421, 9428, 9435, 9442, 9449, 9456, 9463, 9470,
9477, 9484, 9491, 9498, 9505, 9512, 9519, 9526, 9533, 9540, 9547,
9554, 9561, 9568, 9575, 9582, 9589, 9596, 9603, 9610, 9617, 9624,
9631, 9638, 9645, 9652, 9659, 9666, 9673, 9680, 9687, 9694, 9701,
9708, 9715, 9722, 9729, 9736, 9743, 9750, 9757, 9764, 9771, 9778,
9785, 9792, 9799, 9806, 9813, 9820, 9827, 9834, 9841, 9848, 9855,
9862, 9869, 9876, 9883, 9890, 9897, 9904, 9911, 9918, 9925, 9932,
9939, 9946, 9953, 9960, 9967, 9974, 9981, 9988, 9995, 10002,
10009, 10016, 10023, 10030, 10037, 10044, 10051, 10058, 10065,
10072, 10079, 10086, 10093, 10100, 10107, 10114, 10121, 10128,
10135, 10142, 10149, 10156, 10163, 10170, 10177, 10184, 10191,
10198, 10205, 10212, 10219, 10226, 10233, 10240, 10247, 10254,
10261, 10268, 10275, 10282, 10289, 10296, 10303, 10310, 10317,
10324, 10331, 10338, 10345, 10352, 10359, 10366, 10373, 10380,
10387, 10394, 10401, 10408, 10415, 10422, 10429, 10436, 10443,
10450, 10457, 10464, 10471, 10478, 10485, 10492, 10499, 10506,
10513, 10520, 10527, 10534, 10541, 10548, 10555, 10562, 10569,
10576, 10583, 10590, 10597, 10604, 10611, 10618, 10625, 10632,
10639, 10646, 10653, 10660, 10667, 10674, 10681, 10688, 10695,
10702, 10709, 10716, 10723, 10730, 10737, 10744, 10751, 10758,
10765, 10772, 10779, 10786, 10793, 10800, 10807, 10814, 10821,
10828, 10835, 10842, 10849, 10856, 10863, 10870, 10877, 10884,
10891, 10898, 10905, 10912, 10919, 10926, 10933, 10940, 10947,
10954, 10961, 10968, 10975, 10982, 10989, 10996, 11003, 11010,
11017, 11024, 11031, 11038, 11045, 11052, 11059, 11066, 11073,
11080, 11087, 11094, 11101, 11108, 11115, 11122, 11129, 11136,
11143, 11150, 11157, 11164, 11171, 11178, 11185, 11192, 11199,
11206, 11213, 11220, 11227, 11234, 11241, 11248, 11255, 11262,
11269, 11276, 11283, 11290, 11297, 11304, 11311, 11318, 11325,
11332, 11339, 11346, 11353, 11360, 11367, 11374, 11381, 11388,
11395, 11402, 11409, 11416, 11423, 11430, 11437, 11444, 11451,
11458, 11465, 11472, 11479, 11486, 11493, 11500, 11507, 11514,
11521, 11528, 11535, 11542, 11549, 11556, 11563, 11570, 11577,
11584, 11591, 11598, 11605, 11612, 11619, 11626, 11633, 11640,
11647, 11654, 11661, 11668, 11675, 11682, 11689, 11696, 11703,
11710, 11717, 11724, 11731, 11738, 11745, 11752, 11759, 11766,
11773, 11780, 11787, 11794, 11801, 11808, 11815, 11822, 11829,
11836, 11843, 11850, 11857, 11864, 11871, 11878, 11885, 11892,
11899, 11906, 11913, 11920, 11927, 11934, 11941, 11948, 11955,
11962, 11969, 11976, 11983, 11990, 11997, 12004, 12011, 12018,
12025, 12032, 12039, 12046, 12053, 12060, 12067, 12074, 12081,
12088, 12095, 12102, 12109, 12116, 12123, 12130, 12137, 12144,
12151, 12158, 12165, 12172, 12179, 12186, 12193, 12200, 12207,
12214, 12221, 12228, 12235, 12242, 12249, 12256, 12263, 12270,
12277, 12284, 12291, 12298, 12305, 12312, 12319, 12326, 12333,
12340, 12347, 12354, 12361, 12368, 12375, 12382, 12389, 12396,
12403, 12410, 12417, 12424, 12431, 12438, 12445, 12452, 12459,
12466, 12473, 12480, 12487, 12494, 12501, 12508, 12515, 12522,
12529, 12536, 12543, 12550, 12557, 12564, 12571, 12578, 12585,
12592, 12599, 12606, 12613, 12620, 12627, 12634, 12641, 12648,
12655, 12662, 12669, 12676, 12683, 12690, 12697, 12704, 12711,
12718, 12725, 12732, 12739, 12746, 12753, 12760, 12767, 12774,
12781, 12788, 12795, 12802, 12809, 12816, 12823, 12830, 12837,
12844, 12851, 12858, 12865, 12872, 12879, 12886, 12893, 12900,
12907, 12914, 12921, 12928, 12935, 12942, 12949, 12956, 12963,
12970, 12977, 12984, 12991, 12998, 13005, 13012, 13019, 13026,
13033, 13040, 13047, 13054, 13061, 13068, 13075, 13082, 13089,
13096, 13103, 13110, 13117, 13124, 13131, 13138, 13145, 13152,
13159, 13166, 13173, 13180, 13187, 13194, 13201, 13208, 13215,
13222, 13229, 13236, 13243, 13250, 13257, 13264, 13271, 13278,
13285, 13292, 13299, 13306, 13313, 13320, 13327, 13334, 13341,
13348, 13355, 13362, 13369, 13376, 13383, 13390, 13397, 13404,
13411, 13418, 13425, 13432, 13439, 13446, 13453, 13460, 13467,
13474, 13481, 13488, 13495, 13502, 13509, 13516, 13523, 13530,
13537, 13544, 13551, 13558, 13565, 13572, 13579, 13586, 13593,
13600, 13607, 13614, 13621, 13628, 13635, 13642, 13649, 13656,
13663, 13670, 13677, 13684, 13691, 13698, 13705, 13712, 13719,
13726, 13733, 13740, 13747, 13754, 13761, 13768, 13775, 13782,
13789, 13796, 13803, 13810, 13817, 13824, 13831, 13838, 13845,
13852, 13859, 13866, 13873, 13880, 13887, 13894, 13901, 13908,
13915, 13922, 13929, 13936, 13943, 13950, 13957, 13964, 13971,
13978, 13985, 13992, 13999, 14006, 14013, 14020, 14027, 14034,
14041, 14048, 14055, 14062, 14069, 14076, 14083, 14090, 14097,
14104, 14111, 14118, 14125, 14132, 14139, 14146, 14153, 14160,
14167, 14174, 14181, 14188, 14195, 14202, 14209, 14216, 14223,
14230, 14237, 14244, 14251, 14258, 14265, 14272, 14279, 14286,
14293, 14300, 14307, 14314, 14321, 14328, 14335, 14342, 14349,
14356, 14363, 14370, 14377, 14384, 14391, 14398, 14405, 14412,
14419, 14426, 14433, 14440, 14447, 14454, 14461, 14468, 14475,
14482, 14489, 14496, 14503, 14510, 14517, 14524, 14531, 14538,
14545, 14552, 14559, 14566, 14573, 14580, 14587, 14594, 14601,
14608, 14615, 14622, 14629, 14636, 14643, 14650, 14657, 14664,
14671, 14678, 14685, 14692, 14699, 14706, 14713, 14720, 14727,
14734, 14741, 14748, 14755, 14762, 14769, 14776, 14783, 14790,
14797, 14804, 14811, 14818, 14825, 14832, 14839, 14846, 14853,
14860, 14867, 14874, 14881, 14888, 14895, 14902, 14909, 14916,
14923, 14930, 14937, 14944, 14951, 14958, 14965, 14972), class = "Date"),
cases = c(168L, 199L, 214L, 230L, 267L, 373L, 387L, 443L,
579L, 821L, 1229L, 1014L, 831L, 648L, 257L, 203L, 137L, 78L,
82L, 69L, 45L, 51L, 45L, 63L, 55L, 54L, 52L, 27L, 24L, 12L,
10L, 22L, 42L, 32L, 52L, 82L, 95L, 91L, 104L, 143L, 114L,
100L, 83L, 113L, 145L, 175L, 222L, 258L, 384L, 755L, 976L,
879L, 846L, 1004L, 801L, 799L, 680L, 530L, 410L, 302L, 288L,
234L, 269L, 245L, 240L, 176L, 188L, 128L, 96L, 59L, 63L,
44L, 52L, 39L, 50L, 36L, 40L, 48L, 32L, 39L, 28L, 29L, 16L,
20L, 25L, 25L, 48L, 57L, 76L, 117L, 107L, 91L, 90L, 83L,
76L, 86L, 104L, 101L, 116L, 120L, 185L, 290L, 537L, 485L,
561L, 1142L, 1213L, 1235L, 1085L, 1052L, 987L, 918L, 746L,
620L, 396L, 280L, 214L, 148L, 148L, 94L, 107L, 69L, 55L,
69L, 47L, 43L, 49L, 30L, 42L, 51L, 41L, 39L, 40L, 38L, 22L,
37L, 26L, 40L, 56L, 54L, 74L, 99L, 114L, 114L, 120L, 114L,
123L, 131L, 170L, 147L, 163L, 163L, 160L, 158L, 163L, 124L,
115L, 176L, 171L, 214L, 320L, 507L, 902L, 1190L, 1272L, 1282L,
1146L, 896L, 597L, 434L, 216L, 141L, 101L, 86L, 65L, 55L,
35L, 49L, 29L, 55L, 53L, 57L, 34L, 43L, 42L, 13L, 17L, 20L,
27L, 36L, 47L, 64L, 77L, 82L, 82L, 95L, 107L, 96L, 106L,
93L, 114L, 102L, 116L, 128L, 123L, 212L, 203L, 165L, 267L,
550L, 761L, 998L, 1308L, 1613L, 1704L, 1669L, 1296L, 975L,
600L, 337L, 259L, 145L, 91L, 70L, 79L, 63L, 58L, 51L, 53L,
39L, 49L, 33L, 47L, 56L, 32L, 43L, 47L, 19L, 32L, 18L, 34L,
39L, 63L, 57L, 55L, 69L, 76L, 103L, 99L, 108L, 131L, 113L,
106L, 122L, 138L, 136L, 175L, 207L, 324L, 499L, 985L, 1674L,
1753L, 1419L, 1105L, 821L, 466L, 274L, 180L, 143L, 82L, 101L,
72L, 55L, 71L, 50L, 33L, 26L, 25L, 27L, 21L, 24L, 24L, 20L,
18L, 18L, 25L, 23L, 13L, 10L, 16L, 9L, 12L, 16L, 25L, 31L,
36L, 40L, 36L, 47L, 32L, 46L, 75L, 63L, 49L, 90L, 83L, 101L,
78L, 79L, 98L, 131L, 83L, 122L, 179L, 334L, 544L, 656L, 718L,
570L, 323L, 220L, 194L, 125L, 95L, 77L, 46L, 42L, 29L, 35L,
21L, 29L, 16L, 14L, 19L, 15L, 19L, 18L, 21L, 10L, 14L, 7L,
7L, 5L, 9L, 14L, 11L, 18L, 22L, 39L, 36L, 46L, 44L, 37L,
30L, 39L, 37L, 45L, 71L, 59L, 57L, 80L, 68L, 88L, 72L, 74L,
208L, 357L, 621L, 839L, 964L, 835L, 735L, 651L, 400L, 292L,
198L, 85L, 64L, 41L, 40L, 23L, 18L, 14L, 22L, 9L, 19L, 8L,
14L, 12L, 15L, 14L, 4L, 6L, 7L, 7L, 8L, 13L, 10L, 19L, 17L,
20L, 22L, 40L, 37L, 45L, 34L, 26L, 35L, 67L, 49L, 77L, 82L,
80L, 104L, 88L, 49L, 73L, 113L, 142L, 152L, 206L, 293L, 513L,
657L, 919L, 930L, 793L, 603L, 323L, 202L, 112L, 55L, 31L,
27L, 15L, 15L, 6L, 13L, 21L, 10L, 11L, 9L, 8L, 11L, 7L, 5L,
1L, 4L, 7L, 2L, 6L, 12L, 14L, 21L, 29L, 32L, 26L, 22L, 44L,
39L, 47L, 44L, 93L, 145L, 289L, 456L, 685L, 548L, 687L, 773L,
575L, 355L, 248L, 179L, 129L, 122L, 103L, 72L, 72L, 36L,
26L, 31L, 12L, 14L, 14L, 14L, 7L, 8L, 2L, 7L, 8L, 9L, 26L,
10L, 13L, 13L, 5L, 5L, 3L, 6L, 1L, 10L, 6L, 7L, 17L, 12L,
21L, 32L, 29L, 18L, 22L, 24L, 38L, 52L, 53L, 73L, 49L, 52L,
70L, 77L, 95L, 135L, 163L, 303L, 473L, 823L, 1126L, 1052L,
794L, 459L, 314L, 252L, 111L, 55L, 35L, 14L, 30L, 21L, 16L,
9L, 11L, 6L, 6L, 8L, 9L, 9L, 10L, 15L, 15L, 11L, 6L, 3L,
8L, 4L, 7L, 7L, 13L, 10L, 23L, 24L, 36L, 25L, 34L, 37L, 46L,
39L, 37L, 55L, 65L, 54L, 60L, 82L, 55L, 53L, 61L, 52L, 75L,
92L, 121L, 170L, 199L, 231L, 259L, 331L, 357L, 262L, 154L,
77L, 34L, 41L, 21L, 17L, 16L, 7L, 15L, 11L, 7L, 5L, 6L, 13L,
7L, 6L, 8L, 7L, 1L, 11L, 9L, 3L, 9L, 9L, 8L, 15L, 19L, 16L,
10L, 12L, 26L, 35L, 35L, 41L, 34L, 30L, 36L, 43L, 23L, 55L,
107L, 141L, 217L, 381L, 736L, 782L, 663L, 398L, 182L, 137L,
79L, 28L, 26L, 16L, 14L, 8L, 4L, 4L, 6L, 6L, 11L, 4L, 5L,
7L, 7L, 6L, 8L, 2L, 3L, 3L, 1L, 1L, 3L, 3L, 2L, 8L, 8L, 11L,
10L, 11L, 8L, 24L, 25L, 25L, 33L, 36L, 51L, 61L, 74L, 92L,
89L, 123L, 402L, 602L, 524L, 494L, 406L, 344L, 329L, 225L,
136L, 136L, 84L, 55L, 55L, 42L, 19L, 28L, 8L, 7L, 2L, 7L,
6L, 4L, 3L, 5L, 3L, 3L, 0L, 1L, 2L, 3L, 2L, 1L, 2L, 2L, 9L,
4L, 9L, 10L, 18L, 15L, 13L, 12L, 10L, 19L, 15L, 22L, 23L,
34L, 43L, 53L, 47L, 57L, 328L, 552L, 787L, 736L, 578L, 374L,
228L, 161L, 121L, 96L, 58L, 50L, 37L, 14L, 9L, 6L, 15L, 12L,
9L, 1L, 6L, 4L, 7L, 7L, 3L, 6L, 9L, 15L, 22L, 28L, 34L, 62L,
54L, 75L, 65L, 58L, 57L, 60L, 37L, 47L, 60L, 89L, 90L, 193L,
364L, 553L, 543L, 676L, 550L, 403L, 252L, 140L, 125L, 99L,
63L, 63L, 76L, 85L, 68L, 67L, 38L, 25L, 24L, 11L, 9L, 9L,
4L, 8L, 4L, 6L, 5L, 2L, 6L, 4L, 4L, 1L, 5L, 4L, 1L, 2L, 2L,
2L, 2L, 3L, 4L, 4L, 7L, 5L, 2L, 10L, 11L, 17L, 11L, 16L,
15L, 11L, 12L, 21L, 20L, 25L, 46L, 51L, 90L, 123L)), .Names = c("date",
"cases"), row.names = c(NA, -835L), class = "data.frame")