R code prepared for improving graphs & transparency of a colleague's scientometrics analysis, published at https://olivre.com.br/cientometria-da-astronomia [IN PORTUGUESE, by Marcelo Hermes de Lima]

Script_astronomia_MHL2.txt

##############Script R escrito por Eduardo G P Fox em 12/11/2019#########################
#Relativo a um artigo enfatizando o impacto da ciencia nacional em Astrologia, relativo a outros países #relacionados 
#Publicado em: https://olivre.com.br/cientometria-da-astronomia
#idealizado pelos Profs. Marcelo Hermes de Lima (UnB) e Carlos Alberto de Oliveira Torres
#Por favor enviem comentarios sobre este script para meu site pessoal 
####(https://ofoxofox.wixsite.com/research) e/ou Twitter @ofoxofox######

#Pacotes R utilizados
library(ggplot2)
library(stringr)
library(plyr)
library(reshape2)
library(tidyverse)
library(cowplot)
library(magick)
library(TeachingDemos)
#Caso faltar no GUI local, utilizar install.packages() para rodar o script abaixo

#Caso caracteres em Português não estiverem mostrando corretamente, tem que mudar o locale
# no sistema operacional. No meu MacBook OSX 10.12.5 funcionou digitar: 
Sys.setlocale("LC_ALL", "pt_BR.UTF-8")


#Dados recuperados pelo Prof. Marcelo Hermes de Lima a partir da base www.scimagojr.com/ 
#usando a aba "Country Rank", principalmente Citations per Paper (CPP) 

#Inserindo os dados no script R
#Rank Score de países selecionados
#Nota: MHL decidiu por remover deste um valor "outlier": onde diz NA em (AR) era "8.4"
Rank_score<-data.frame(
Anos=c(2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018),
BR=c(1.5, 0.9, 1.7, 3.2, 2.4, 3.6, 1.3, 3.95, 5, 6.05, 5),
CH= c(5.6, 5.7, 5.6, 5.9, 8.4, 6.1, 4.9, 5.8, 4.7, 7.1, 6.75),
PT= c(NA, 6.3, 3.9, 6.8, 5.7, 5, 6.7, 4.7, 5.3, 7.6, 7.75),
AR= c(1.1, 1.4, 2.5, 1.6, 3.5, 1.9, NA, 1.8, 5.7, 3.7, 1.5)
)

#% da média dos 3 países líderes na área nestes anos
Perc_top_3<-data.frame(
Anos=c(2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018),
BR=c(40, 31.7, 37.7, 53.5, 53.1, 60.6, 27, 72.4, 49, 64.9, 57.2),
CH=c(64, 71.7, 76.1, 74.7, 90.8, 80.2, 47.8, 80, 46.8, 68.5, 64.6),
PT=c(NA, NA, NA, NA, NA, NA, NA, 75.3, 43.6, 71.3, 75.6),
AR=c(NA, NA, NA, NA, NA, NA, 170.3, 43.5, 43.8, 55.3, 27.7)
)

#dados somente do Brasil
BR_CPP_Rel<-data.frame(
Anos=c(2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018),
top_3=c(53.5, 53.1, 60.6, 27, 72.4, 49, 64.9, 57.2),
mundo=c(77.8, 72.7, 83.4, 56.7, 98.3, 114.9, 114.8, 107.5)
)

#Atenção: MHL retirou um ponto outlier abaixo, logo NA (Impacto) era "27.0" no original
Brasil<-data.frame(
Anos=seq(2007, 2019, 1),
Quantidade=c(NA,327, 411, 371, 411, 399, 416, 560, 559, 643, 568, 669, NA),
Impacto=c(NA,40.0, 31.7, 37.7, 53.5, 53.1, 60.6, NA, 72.4, 49.0, 64.9, 57.2, NA)
)

CPP<-data.frame(
Anos=c(2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018),
mundo=c(31.21, 26.73, 25.57, 24.5, 15.99, 16.47, 7.2, 2.075),
top_3=c(45.4, 36.56, 35.17, 51.7, 21.7, 38.6, 12.75, 3.9)
)

R_Score<-data.frame(
Anos=seq(2000, 2020, 1),
BR=c(rep(NA, 3), 1.2, 2, 1.5, 2.7, 3.2, 0.9, 1.5, 1.5, 0.9, 1.7, 3.2, 2.4, 3.6, 1.3, 3.95, 5, 6.05, 5),
CH=c(rep(NA, 3), 6.1, 8, 9.3, 6.3, 4.2, 5.6, 6.5, 5.6, 5.7, 5.6, 5.9, 8.4, 6.1, 4.9, 5.8, 4.7, 7.1, 6.75)
)

Quant_Art<-data.frame(
Anos=seq(2007, 2019, 1),
BR=c(NA, 327, 411, 371, 411, 399, 416, 560, 559, 643, 568, 669, NA),
CH=c(NA, 503, 561, 548, 705, 688, 781, 964, 921, 1023, 967, 1060, NA),
AR=c(NA, 206, 191, 212, 257, 188, 197, 202, 196, 236, 208, 217, NA)
)


#Plotando as Figuras como no original, segundo a estética pedida pelo primeiro autor MHL

#Imagens abaixo obtidas a partir do website https://www.flaticon.com/packs/countrys-flags
#Icons made by <a href="https://www.flaticon.com/authors/freepik" title="Freepik">Freepik</a> from <a href="https://www.flaticon.com/" title="Flaticon">www.flaticon.com</a></div>
#Ajustando ícones de bandeiras
brazil_icon<-image_read('/Users/egoncal2/Downloads/brazil.png')
chile_icon<-image_read('/Users/egoncal2/Downloads/chile.png')
argentina_icon<-image_read('/Users/egoncal2/Downloads/argentina.png')
portugal_icon<-image_read('/Users/egoncal2/Downloads/portugal.png')

#Fig. 1 mostra nuvens de pontos (obtidas por erro padrão) e tendências de diferentes países
#Versão da Figura 1A usando regressão em curva (em método "default") 
Rank1c<-ggplot(R_Score)+
geom_point(aes(y=BR, x=Anos), colour="#009C3B", size = 3)+
geom_smooth(aes(y=BR, x=Anos), size = 1, colour="#009C3B", linetype = "dashed") +
geom_point(aes(y=CH, x=Anos), colour="#0039A6", size = 3)+
geom_smooth(aes(y=CH, x=Anos), size = 1, colour="#0039A6", linetype = "dashed") +
scale_y_continuous(limits=c(0,10), breaks=c(2,4,6,8,10), expand = c(0, 0))+
labs(title = "Rank Score", x = "Ano de Publicação", y="Rank-Score")+
theme_cowplot(12)+
theme(
legend.position = "bottom",
plot.title = element_text(color="red", size=16, face="bold.italic"),
axis.title.x = element_text(size=14, face="bold"),
axis.title.y = element_text(size=14, face="bold",color="white"),
axis.text.x=element_text(size=14),
axis.text.y=element_text(size=12))

#Ajustando icones de imagem
Rank1ic<-ggdraw(Rank1c) +
draw_image(brazil_icon, scale=0.1, y=-0.22, x=-0.3)+
draw_image(chile_icon, scale=0.1, y=0.22, x=-0.3)

#Figura 1B do artigo aplicando uma regressão linear simples
Rank2<-ggplot(Rank_score[,c(1,4:5)])+
geom_point(aes(y=PT, x=Anos), colour="#ff0000", size = 3)+
geom_smooth(aes(y=PT, x=Anos), method = lm, size = 1, span = 1, colour="#ff0000", linetype = "dashed") +
geom_point(aes(y=AR, x=Anos), colour="#75aadb", size = 3)+
geom_smooth(aes(y=AR, x=Anos), method = lm, size = 1, span = 1, colour="#75aadb", linetype = "dashed") +
scale_x_continuous(limits=c(2006,2019), breaks=c(2009,2011,2013,2015,2017), expand = c(0, 0))+
scale_y_continuous(limits=c(0,10), breaks=c(2,4,6,8,10), expand = c(0, 0))+
theme_cowplot(12)+
labs(title = " ", x = "Ano de Publicação", y="Rank-Score")+
theme(
legend.position = "bottom",
plot.title = element_text(color="red", size=16, face="bold.italic"),
axis.title.x = element_text(size=14, face="bold"),
axis.title.y = element_text(size=14, face="bold",color="white"),
axis.text.x=element_text(size=14),
axis.text.y=element_text(size=12))

#Ajustando ícones de imagem
Rank2i<-ggdraw(Rank2) +
draw_image(portugal_icon, scale=0.1, x=-0.28, y=0.05)+
draw_image(argentina_icon, scale=0.1, y=-0.2, x=-0.28)

#Plotando o conjunto sem detalhes 
plot_grid(Rank1ic, Rank2i)

#Figura 2 segundo estética pedida pelo Prof. MHL
Fig2<-ggplot(Perc_top_3[,1:3])+
geom_point(aes(y=BR, x=Anos), colour="#009C3B", size = 3)+
geom_smooth(aes(y=BR, x=Anos), method = "loess", se = FALSE, size = 1, span = 1, colour="#009C3B") +
geom_point(aes(y=CH, x=Anos), colour="#0039A6", size = 3)+
geom_smooth(aes(y=CH, x=Anos), method = "loess", se = FALSE, size = 1, span = 1, colour="#0039A6", linetype = "dotted") +
scale_x_continuous(limits=c(2007,2019), breaks=c(2009,2011,2013,2015,2017), expand = c(0, 0))+
scale_y_continuous(limits=c(0,100), breaks=c(20,40,60,80,100), expand = c(0, 0))+
theme_light()+
labs(title = "% dos Top 3 em CPP", x = "Ano de Publicação")+
theme(
legend.position = "bottom",
plot.title = element_text(color="red", size=16, face="bold.italic"),
axis.title.x = element_text(size=14, face="bold"),
axis.title.y = element_text(size=14, face="bold",color="white"),
axis.text.x=element_text(size=14),
axis.text.y=element_text(size=12))
#Adicionando ícones de países
ggdraw(Fig2) +
draw_image(chile_icon, scale=0.1, x=-0.3, y=0.05)+
draw_image(brazil_icon, scale=0.1, y=-0.2, x=-0.3)

#Figura 3: Quantidade de Artigos por Ano de Diferentes países 
Quant1<-ggplot(Quant_Art)+
geom_point(aes(y=BR, x=Anos), colour="#009C3B", size = 3)+
geom_smooth(aes(y=BR, x=Anos), method = "loess", se = FALSE, size = 1, span = 0.2, colour="#009C3B") +
geom_point(aes(y=CH, x=Anos), colour="#0039A6", size = 3)+
geom_smooth(aes(y=CH, x=Anos), method = "loess", se = FALSE, size = 1, span = 0.2, colour="#0039A6") +
geom_point(aes(y=AR, x=Anos), colour="#75aadb", size = 3)+
geom_smooth(aes(y=AR, x=Anos), method = "loess", se = FALSE, size = 1, span = 0.2, colour="#75aadb") +
scale_x_continuous(limits=c(2006.5,2018.5), breaks=c(2009,2011,2013,2015,2017), expand = c(0, 0))+
theme_light()+
labs(title = "Número de Artigos", x = "Ano de Publicação")+
theme(
legend.position = "bottom",
plot.title = element_text(color="red", size=16, face="bold.italic"),
axis.title.x = element_text(size=14, face="bold"),
axis.title.y = element_text(size=14, face="bold",color="white"),
axis.text.x=element_text(size=14),
axis.text.y=element_text(size=12))
#Adicionando os icones de países
ggdraw(Quant1) +
draw_image(chile_icon, scale=0.09, y=-0.08, x=-0.33)+
draw_image(brazil_icon, scale=0.09, y=-0.21, x=-0.33)+
draw_image(argentina_icon, scale=0.09, y=-0.33, x=-0.33)

#Figura 4 aqui apelando para R "básico", pois ggplot2 dificulta plotar eixos paralelos
#Regressão linear simples para comparar retas obtidas em termos de tendência
Brasil.Impacto<-Brasil$Impacto
Brasil.Quantidade<-Brasil$Quantidade
Brasil.Anos<-Brasil$Anos
Quantidade.lm<-lm(Brasil.Quantidade ~ Brasil.Anos)
Impacto.lm<-lm(Brasil$Impacto ~ Brasil.Anos)

#Para restringir a regressão para determinados pontos foi preciso usar a função clipplot()
par(mar=c(5, 4, 4, 6) + 0.1)
plot(Brasil$Quantidade ~ Brasil$Anos, col = "red", ylim=c(0,1000), xlim=c(2007,2019.5), pch=19, 
cex =2, main="Relacionando Quantidade vs. Impacto", xlab="Anos", ylab=" ", cex.lab=2, cex.main=1.5) 
clipplot(abline(Quantidade.lm, col = "red"), xlim=c(2008,2018))
par(new = T)
plot(Brasil$Impacto ~ Brasil$Anos, col = "blue", ylim=c(0,100), xlim=c(2007,2019.5), 
axes = FALSE, bty = "n", xlab = "", pch=18, cex =2, ylab = "") 
axis(side=4)
clipplot(abline(Impacto.lm, col = "blue"), xlim=c(2008,2018))
mtext("Percentual dos top-3" ,side=4,col="blue",line=2, cex=1.1) 
mtext("Número de Artigos",side=2,col="red",line=2, cex=1.1)
rasterImage(brazil_icon, 2007, 80, 2008.5,100)
legend("bottomright",legend=c("Quantidade","Impacto"),
  text.col=c("red","blue"),pch=c(19,18),col=c("red","blue"))


#### Por favor, alguém poderia sugerir uma implementação a Fig. 4 com estética ggplot2 ? #### 

############################### Obrigado ########################################