Final

Source:Getty images/istockphoto

Introduction:

The dataset originates from the esteemed Pew Research Center and their data is featured in their comprehensive analysis of religious restrictions specifically in the Asian-Pacific, Middle East, and African regions. I had to take their data and compile it into an excel. It includes key variables including country, year, location, Government Restriction Index (GRI), and Social Hostilities Index (SHI).

The Government Restriction Index (GRI) was created through an examination of governmental regulations, policies, and laws that impose limitations on religious beliefs and practices. This includes assessing whether governments provide preferential treatment to specific religious groups, or enact bureaucratic hurdles hindering certain groups from accessing benefits. Based on these factors each country was given a score; the lower the better.

Conversely, the Social Hostilities Index (SHI) was crafted by tracking various events such as incidents of religion-related harassment, mob violence, terrorism or militant activities, and conflicts arising from religious conversions or the display of religious symbols and attire. I am curious to find out the relationship between GRI and SHI along with which countries in these regions have the lowest and highest values for each.

The data set presents an invaluable opportunity to delve into historical trends and discern the intricate relationship between GRI and SHI, enabling insightful cross-continental comparisons to uncover nuances in how religious restrictions impact social hostilities across different regions.

install.packages("highcharter")

Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.4'
(as 'lib' is unspecified)

install.packages("leaflet")

Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.4'
(as 'lib' is unspecified)

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.5.1     ✔ tibble    3.2.1
✔ lubridate 1.9.3     ✔ tidyr     1.3.1
✔ purrr     1.0.2     

── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(readr)
library(ggplot2)
library(highcharter)

Registered S3 method overwritten by 'quantmod':
  method            from
  as.zoo.data.frame zoo 
Highcharts (www.highcharts.com) is a Highsoft software product which is
not free for commercial and Governmental use

library(dplyr)
library(purrr)
library(RColorBrewer)
library(leaflet)
DF <- read_csv("Final Dataset(Sheet1).csv")

Rows: 339 Columns: 7
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (4): Country, Location, GRI, SRI
dbl (3): Year, Lat, Lon

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

DF_filtered <- DF %>% filter(Year == 2016)

DF_filtered$GRI <- as.numeric(DF_filtered$GRI)
DF_filtered$SRI <- as.numeric(DF_filtered$SRI)

DF_filtered_chart <- ggplot(DF_filtered, aes(x = GRI, y = SRI)) +
  xlab("Government Restriction Index ") +
  ylab("Social Hostilities Index") +
  ggtitle("GRI VS SRI In 2016") +
  theme_minimal(base_size = 14, base_family = "serif") +
  scale_x_continuous(breaks = 1:10) +  # Specify the breaks for the x-axis
  scale_y_continuous(breaks = 1:10)    # Specify the breaks for the y-axis

DF_filtered_chart +
 geom_point()

p4 <- DF_filtered_chart + geom_point() + geom_smooth(color = "red")
p4

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

p5 <- DF_filtered_chart + geom_point() + geom_smooth(method='lm',formula=y~x)
p5

model <- lm(GRI ~ SRI, data = DF_filtered)
summary(model)

Call:
lm(formula = GRI ~ SRI, data = DF_filtered)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.9313 -1.8475 -0.6447  1.7849  5.3680 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  2.96197    0.33615   8.811 2.32e-14 ***
SRI          0.36155    0.08747   4.133 7.08e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.331 on 108 degrees of freedom
Multiple R-squared:  0.1366,    Adjusted R-squared:  0.1286 
F-statistic: 17.09 on 1 and 108 DF,  p-value: 7.082e-05

y= 2.96197 + 0.36155 x GRI

P value= 7.082

R^2= 0.1366

The P-value suggests there is no correlation between the SHI and the GRI. R^2 suggests a relatively weak relationship between the two variables.

DF_filtered_chart+
geom_point(size = 3, alpha = 0.5, aes(color = Location)) +
 geom_smooth(method = lm, se =FALSE, color = "black", lty = 2, linewidth = 0.3)

`geom_smooth()` using formula = 'y ~ x'

DF_filtered2 <- DF %>%
  filter(!is.na(Year) & !is.na(Location) & !is.na(GRI) & !is.na(SRI)) %>%
  filter(Year %in% c(2007, 2016))

DF_filtered2$GRI <- as.numeric(DF_filtered2$GRI)

Warning: NAs introduced by coercion

DF_filtered2$SRI <- as.numeric(DF_filtered2$SRI)

Warning: NAs introduced by coercion

DF_filtered_chart2 <- ggplot(DF_filtered2, aes(x = GRI, y = SRI, col = Location)) +
  geom_point(alpha = 0.8) +
  guides(size = "none") +
  theme(plot.title = element_blank(), legend.title = element_blank()) +
  xlab("Government Restriction Index") +
  ylab("Social Hostilities Index") +
  geom_text(aes(x = 6, y = 9, label = Year), cex = 12, color = "grey") +
  facet_grid(. ~ Year) +
  theme(strip.background = element_blank(),
        strip.text.x = element_blank(),
        strip.text.y = element_blank(),
        legend.position = "top") +
  scale_x_continuous(breaks = 0:10, limits = c(0, 10)) +  # Specify the breaks and limits for the x-axis
  scale_y_continuous(breaks = 0:10, limits = c(0, 10))    # Specify the breaks and limits for the y-axis

print(DF_filtered_chart2)

Warning: Removed 1 row containing missing values or values outside the scale range
(`geom_point()`).

hc <- highchart()

unique_locations <- na.omit(unique(DF_filtered$Location))  # Remove NA values from unique locations

# Define colors for each location
location_colors <- c("#00FF2E", "#057DFF", "yellow")  

# Add series for each location, excluding Series 4
for (i in seq_along(unique_locations)) {
  loc <- unique_locations[i]
  
  # Skip if loc is NA
  if (is.na(loc)) {
    next  # Skip to the next iteration
  }
  
  data <- DF_filtered[DF_filtered$Location == loc, ]
  
  hc <- hc %>%
    hc_add_series(
      data = data,
      type = "scatter",
      hcaes(x = GRI, y = SRI, text = paste("Location:", Location, "<br>Country:", Country, "<br>SRI:", SRI, "<br>GRI:", GRI)),
      name = loc,
      color = location_colors[i],  # Assign color based on location
      marker = list(symbol = "circle", radius = 4)
    )
}

hc <- hc %>%
  hc_tooltip(pointFormat = "<b>{point.text}</b>") %>%
  hc_xAxis(title = list(text = "Government Restriction Index")) %>%  # Add x-axis label
  hc_yAxis(title = list(text = "Social Hostilities Index")) %>%  # Add y-axis label
  hc_title(text = "Scatter Plot of GRI vs. SHI In 2016") %>%  # Add chart title
  hc_add_theme(hc_theme_darkunica()) %>%  # Add Dark Unica theme
  hc_chart(events = list(load = JS("function() {
    this.renderer.text('Source: Pew Research Center', this.plotLeft + 10, this.plotTop + this.plotHeight + 40)
      .attr({ zIndex: 6, fill: 'white' }).add();
  }")))
  
print(hc) 

This scatterplot shows the GRI and SHI for countries in the Asia-Pacific, Middle East- North Africa, and Sub-Saharan Africa. The tooltip shows the country and the specific SHI and GRI for each point for the year 2016.

Top10df <- DF_filtered %>%
  filter(GRI >= 6.7 & GRI <= 10) %>%
  top_n(10, GRI)

Top10df <- as.data.frame(Top10df)  

Top10df_reshaped <- Top10df %>%
  select(Country, GRI, SRI) %>%
  pivot_longer(cols = c(GRI, SRI), names_to = "Index_Type", values_to = "Value")

hc <- highchart() %>% 
  hc_xAxis(categories = unique(Top10df_reshaped$Country),
           title = list(text = "Country")) %>% 
  hc_add_series(name = "GRI", 
                data = Top10df_reshaped %>% filter(Index_Type == "GRI") %>% pull(Value)) %>% 
  hc_add_series(name = "SRI", 
                data = Top10df_reshaped %>% filter(Index_Type == "SRI") %>% pull(Value)) %>%
  hc_title(text = "Top 10 Countries with the Highest Government Restriction Index (GRI) in 2016") %>%
  hc_yAxis(title = list(text = "Value"))

hc

hc <- hc %>% 
  hc_chart(borderColor = '#EBBA95',
           borderRadius = 10,
           borderWidth = 2,
           backgroundColor = list(
             linearGradient = c(0, 0, 500, 500),
             stops = list(
               list(0, 'rgb(255, 255, 255)'),
               list(1, 'rgb(200, 200, 255)')
             ))) %>%
  hc_xAxis(title = list(text = "Country")) %>%
  hc_yAxis(title = list(text = "Value")) %>%
  hc_title(text = "Top 10 Countries with the Highest Government Restriction Index (GRI) in 2016")

print(hc)

# Reshape the data to have separate columns for each country's GRI and SRI values
Top10df_reshaped <- Top10df |>
  select(Country, GRI, SRI,Location) |>
  pivot_longer(cols = c(GRI, SRI), names_to = "Index_Type", values_to = "Value")

# Create Highcharter plot with 3D columns
hc <- highchart() |> 
  hc_chart(type = "column",
           options3d = list(enabled = TRUE, beta = 15, alpha = 15)) |>
  hc_xAxis(categories = unique(Top10df_reshaped$Country),
           title = list(text = "Country", margin = 20)) |>  # Adjust margin for the x-axis title
  hc_add_series(name = "GRI", 
                data = Top10df_reshaped |> filter(Index_Type == "GRI") |> pull(Value)) |> 
  hc_add_series(name = "SRI", 
                data = Top10df_reshaped |> filter(Index_Type == "SRI") |> pull(Value)) |>
  hc_title(text = "GRI and SRI for Top 10 Countries") |>
  hc_yAxis(title = list(text = "Value"))
hc

color_ranges <- list(
  list(from = 0, to = 1.4, color = "#90EE90"),  # Light green for values less than or equal to 1.4
  list(from = 1.5, to = 3.5, color = "#FFFF00"),   # Yellow for values between 1.5 and 3.5
  list(from = 3.6, to = 7.1, color = "#FFA500"),   # Orange for values between 3.6 and 7.1
  list(from = 7.2, to = 10, color = "#8B0000")    # Dark red for values 7.2 or higher
)

gri_color <- "#A9A9A9"  # Light grey color for GRI series
sri_color <- "#454545"  # Blue color for SRI series

hc <- highchart() |>
  hc_chart(type = "column") |>
  hc_xAxis(categories = unique(Top10df_reshaped$Country),
           title = list(text = "Country", margin = 20)) |>
  hc_add_series(name = "GRI", 
                data = Top10df_reshaped |>
                  filter(Index_Type == "GRI") |>
                  pull(Value),
                color = gri_color) |> # Set color for GRI series
  hc_add_series(name = "SHI", 
                data = Top10df_reshaped |>
                  filter(Index_Type == "SRI") |>
                  pull(Value),
                color = sri_color) |> # Set color for SRI series
  hc_title(text = "Top 10 Countries with the Highest Government Restriction Index (GRI) in 2016") |>
  hc_yAxis(title = list(text = "Value"),  # Configure y-axis title
           plotBands = color_ranges) 

hc <- hc |> hc_chart(events = list(load = JS("function() {
    this.renderer.text('Source: Pew Research Center', this.plotLeft + 10, this.plotTop + this.plotHeight + 100)
      .attr({ zIndex: 6, fill: 'black' }).add();
  }")))

hc

hc <- highchart() |>
  hc_chart(type = "column",
           options3d = list(enabled = TRUE, beta = 15, alpha = 15)) %>%
  hc_xAxis(categories = unique(Top10df_reshaped$Country),
           title = list(text = "Country", margin = 20)) %>%
  hc_add_series(name = "GRI", 
                data = Top10df_reshaped %>%
                  filter(Index_Type == "GRI") %>%
                  pull(Value),
                color = gri_color) %>%  # Set color for GRI series
  hc_add_series(name = "SHI", 
                data = Top10df_reshaped %>%
                  filter(Index_Type == "SRI") %>%
                  pull(Value),
                color = sri_color) %>% # Set color for SRI series
  hc_title(text = "Top 10 Countries with the Highest Government Restriction Index (GRI) in 2016") %>%
  hc_yAxis(title = list(text = "Value"),  # Configure y-axis title
         plotBands = color_ranges) 

hc <- hc %>% hc_chart(events = list(load = JS("function() {
    this.renderer.text('Source: Pew Research Center', this.plotLeft + 10, this.plotTop + this.plotHeight + 40)
      .attr({ zIndex: 6, fill: 'black' }).add();
  }")))

hc

This bar graph shows the values of GRI and SHI in the top 10 countries with the highest Government Restriction Index. It also shows the bar graph for SHI. The background colors represent how Pew Research Center defines low, mild, high, and very high restrictions in regards to SHI. This way we can visually compare the heights of the SHI bars to see how they rank according to Pew Research Center.

While exploring the Pew Research Center website they discussed how Maldives, where Islam is the state religion, non-muslims groups are forbidden to make houses of worship. Similarly, under Egyptian law only allows membersof recognized Religous groups (Sunni Islam, Christianity, and Judaism) to create places of worship and to express their faith in public.

top_10_highest_sri_countries <- DF_filtered |>
  arrange(desc(SRI)) |>  # Arrange the dataframe by SRI in descending order
  head(10)  # Get the top 100 rows <- DF_filtered %>%

leaflet(top_10_highest_sri_countries) %>%
  addTiles() %>%  # Add default OpenStreetMap tiles
  addMarkers(~Lon, ~Lat, popup = ~Country)

palette <- colorNumeric(palette = "RdYlBu", domain = top_10_highest_sri_countries$SRI)

leaflet(top_10_highest_sri_countries) %>%
  addTiles() %>%  # Add default OpenStreetMap tiles
  addCircleMarkers(~Lon, ~Lat, 
                   radius = 8,  # Increased radius for larger points
                   color = ~palette(SRI),
                   fillOpacity = 0.7,
                   popup = ~paste("SRI: ", SRI)) %>%
  addLegend("bottomright", pal = palette, values = ~SRI, title = "Value", opacity = 1) %>%
  addControl(
    position = "topleft",
    html = "<h3>Top 10 Countries With The Highest SHI values in 2016</h3>"
  ) %>%
  addControl(
    position = "bottomleft",
    html = "<p>Pew Research Center</p>"
  )

Finally, the last chart shows a map and the ranking of the top 10 Countries with the highest SHI in 2016 in the regions I specified earlier. Out of the top ten the value gradiant shows how the countries rank among themselves with red being the lowest and dark blue being the highest.