Primero, un dput de mi base de datos:
> dput(antib)
structure(list(X = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13,
14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29,
30, 31, 32, 33, 34, 35, 36, 37, 38), state = structure(c(2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L), .Label = c("Control", "WV"), class = "factor"),
year = c(2001, 2001, 2001, 2002, 2002, 2002, 2002, 2003,
2003, 2003, 2003, 2004, 2004, 2004, 2004, 2005, 2005, 2005,
2005, 2001, 2001, 2001, 2002, 2002, 2002, 2002, 2003, 2003,
2003, 2003, 2004, 2004, 2004, 2004, 2005, 2005, 2005, 2005
), quarter = c(2, 3, 4, 1, 2, 3, 4, 1, 2, 3, 4, 1, 2, 3,
4, 1, 2, 3, 4, 2, 3, 4, 1, 2, 3, 4, 1, 2, 3, 4, 1, 2, 3,
4, 1, 2, 3, 4), market_share = c(0.52714, 0.49988, 0.49205,
0.49713, 0.49724, 0.49729, 0.49972, 0.50085, 0.46391, 0.4343,
0.42327, 0.41887, 0.38585, 0.36963, 0.3565, 0.33753, 0.31054,
0.30456, 0.28855, 0.48893, 0.47956, 0.47722, 0.48834, 0.47811,
0.47284, 0.47466, 0.47867, 0.47851, 0.48252, 0.49165, 0.49599,
0.48444, 0.47741, 0.47537, 0.46563, 0.45466, 0.45886, 0.45484
), groups = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("Tratamiento",
"Control"), class = "factor"), time = c(1, 2, 3, 4, 5, 6,
7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 1, 2, 3,
4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19),
level = c(0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1), trend = c(0, 0, 0, 0, 0, 0, 0, 0, 1, 2, 3, 4, 5,
6, 7, 8, 9, 10, 11, 0, 0, 0, 0, 0, 0, 0, 0, 1, 2, 3, 4, 5,
6, 7, 8, 9, 10, 11), trat_time = c(1, 2, 3, 4, 5, 6, 7, 8,
9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), trat_level = c(0,
0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), trat_trend = c(0,
0, 0, 0, 0, 0, 0, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), predicted = structure(c(`1` = 0.510175716539688,
`2` = 0.507834785635212, `3` = 0.505493854730736, `4` = 0.50315292382626,
`5` = 0.500811992921784, `6` = 0.498471062017308, `7` = 0.496130131112831,
`8` = 0.493789200208355, `9` = 0.459358734732496, `10` = 0.441916641822141,
`11` = 0.424474548911787, `12` = 0.407032456001432, `13` = 0.389590363091078,
`14` = 0.372148270180723, `15` = 0.354706177270369, `16` = 0.337264084360015,
`17` = 0.31982199144966, `18` = 0.302379898539306, `19` = 0.284937805628951,
`20` = 0.484548114140817, `21` = 0.483373084731565, `22` = 0.482198055322314,
`23` = 0.481023025913063, `24` = 0.479847996503811, `25` = 0.47867296709456,
`26` = 0.477497937685308, `27` = 0.476322908276057, `28` = 0.491051942381428,
`29` = 0.48766325893013, `30` = 0.484274575478832, `31` = 0.480885892027533,
`32` = 0.477497208576235, `33` = 0.474108525124936, `34` = 0.470719841673638,
`35` = 0.46733115822234, `36` = 0.463942474771041, `37` = 0.460553791319743,
`38` = 0.457165107868445), label = "Fitted values"), contrafactual = c(0,
0, 0, 0, 0, 0, 0, 0, 0.507352332818502, 0.502797747871979,
0.498243162925456, 0.493688577978932, 0.489133993032409,
0.484579408085886, 0.480024823139363, 0.47547023819284, 0.470915653246317,
0.466361068299793, 0.46180648335327, 0, 0, 0, 0, 0, 0, 0,
0, 0.475147878866806, 0.473972849457554, 0.472797820048303,
0.471622790639051, 0.4704477612298, 0.469272731820549, 0.468097702411297,
0.466922673002046, 0.465747643592794, 0.464572614183543,
0.463397584774292)), class = "data.frame", row.names = c(NA,
-38L))
Me interesa que mi gráfico tenga 6 elementos:
1. Valores observados para el grupo de tratamiento (GT)
2. Valores observados para el grupo de control (GC)
3. Valores predichos para el GT
4. Valores predichos para el GC
5. Contrafactual para el GT
6. Contrafactual para el GC
Los primeros dos están en la columna market_share, los predichos en "predicted" y los últimos en "contrafactual".
Ahora bien, ya tengo un avance creo que importante:
ggplot(antib, aes(x=paste0(year, "-", quarter),
y=market_share, group= groups, color= groups)) +
geom_line(size=1.2) +
labs(x = "", y="", title="Cuota de mercado para grupos de tratamiento y control") +
theme(legend.title=element_blank()) +
geom_vline(xintercept = 8.5, size=1) +
geom_point(aes(y = antib$predicted, col = groups)) +
geom_segment(aes(x = antib$time[1], y = antib$predicted[1], xend = antib$time[8], yend = antib$predicted[8]))
El cual nos da la siguiente imagen:
Sin embargo, los valores predichos para el grupo de tratamiento están conectados entre sí por una línea que no corresponde al color asignado para dicho grupo y no entiendo por qué.
Si intento:
ggplot(antib, aes(x=paste0(year, "-", quarter),
y=market_share, group= groups, color= groups)) +
geom_line(size=1.2) +
labs(x = "", y="", title="Cuota de mercado para grupos de tratamiento y control") +
theme(legend.title=element_blank()) +
geom_vline(xintercept = 8.5, size=1) +
geom_point(aes(y = antib$predicted, col = groups)) +
geom_segment(aes(x = antib$time[1], y = antib$predicted[1], xend = antib$time[8], yend = antib$predicted[8],
colour=antib$groups[antib$groups=="Tratamiento"]))
Especificando como color para el segmento el color asignado al grupo de tratamiento, obtengo este error:
Error: Aesthetics must be either length 1 or the same as the data (38): x, y, xend, yend, colour, group
Debo señalar que me interesa que los valores predichos aparezcan por segmentos para los períodos pre-post intervención (señalados con la línea vertical negra) y que, en general, todos los segmentos de línea tengan el mismo color según el grupo del que se trate.
En otras palabras, la gráfica que busco con ggplot es similar a esta:
En la que, como puede apreciarse, las líneas y segmentos tienen colores más o menos similares para distinguir entre GT y GC.
Cualquier ayuda u orientación será muy apreciada.