library(tidyverse)
library(ggmagnify)

pca <- read.table("https://raw.githubusercontent.com/paulobarros/assets/main/pcs.txt", col.names = c("ind","PC1","PC2","PC3","gg")) |>
    select(-ind) |>
    mutate(gg = factor(gg)) |>
    mutate(colors = case_when(
        gg == "BEL" ~ alpha("#124559",.75),
        gg == "BRO" ~ alpha("#e0aaff",.75),
        gg == "CC" ~ alpha("#00b4d8",.75),
        gg == "China" ~ alpha("#80ed99",.75),
        gg == "WEL" ~ alpha("#bb9457",.75),
        gg == "TT" ~ alpha("#ff006e",.75),
        gg == "Nigeria" ~ alpha("#fca311",.75),
        gg == "Indio" ~ alpha("#6c757d",.75),
        gg == "GG" ~ alpha("#f07167",.75)
    ),
    lines = case_when(
        gg == "BEL" ~ alpha("#124559",.95),
        gg == "BRO" ~ alpha("#e0aaff",.95),
        gg == "CC" ~ alpha("#00b4d8",.95),
        gg == "China" ~ alpha("#80ed99",.95),
        gg == "WEL" ~ alpha("#bb9457",.95),
        gg == "TT" ~ alpha("#ff006e",.95),
        gg == "Nigeria" ~ alpha("#fca311",.95),
        gg == "Indio" ~ alpha("#6c757d",.95),
        gg == "GG" ~ alpha("#f07167",.95)
    ))

ggColors <- pca$colors |> unique()
ggLines <- pca$lines |> unique()

from <- c(xmin = -0.015,
          xmax = 0.015,
          ymin = -0.065,
          ymax = -0.02)


to <- c(0.05, 0.2, -0.04,0.03)

inset <- pca |> 
filter(!gg %in% c("TT","China")) |>
ggplot(aes(x = PC1, y = PC2,)) +
geom_point(aes(fill = gg),
color = alpha("black",.75),
shape = 21,
size = 6,
show.legend = FALSE) +
scale_fill_manual(values = ggColors |> unique()) +
labs(x = element_blank(), y  = element_blank()) +
#scale_color_manual(values = ggLines |> unique()) +
theme_bw(base_family = "Barlow Condensed", base_size = 16) +
scale_x_continuous(breaks = (seq(-0.005, 0.001, 0.002)), limits = c(-0.005,0.0011)) +
scale_y_continuous(breaks = (seq(-0.06, -0.023, 0.015)), limits = c(-0.06,-0.023))

pca |> 
#filter(!gg %in% c("TT","China")) |>
ggplot(aes(x = PC1, y = PC2,)) +
geom_point(aes(fill = gg),
color = alpha("black",.75),
shape = 21,
size = 6) +
scale_fill_manual(values = ggColors |> unique()) +
labs(x = "PC1 (10.34%)",
     y = "PC2 (7.73%)",) +
#scale_color_manual(values = ggLines |> unique()) +
theme_bw(base_family = "Barlow Condensed", base_size = 16) +
geom_magnify(from = from,
            to = to,
            shadow = TRUE,
            plot = inset,
            axes = "xy"
            )

library(tidyverse)
library(cowplot)
library(ggpubr)
library(ggmagnify)
library(emmeans)

# Simulate data ----
n <- 1000                    # don't change this unless necessary (plots might be fragile)
set.seed(1234)
x  <-  rnorm(n, 486, 142)    # generate macular hole inner opening data with mean 486 and sd = 152
z  <-  10.89 - 0.016 * x     # generate variable that is linear combination of intercept = 10.89 and coefficient for macular hole -0.016 (logit scale)
pr  <-  1/(1 + exp(-z))      # generate probabilities from this
y  <-  rbinom(n, 1, pr)      # generate outcome variable as a function of those probabilities

# Create dataframe from these:
df <-  data.frame(y = y, x = x, z = z, pr = pr)
df <- df |> 
  filter(x > 100) # only include those with thickness > 100

# Logistic regression model ----
summary(mod_logistic <- glm(y ~ x, data = df, family = "binomial"))

# Predictions ----
# Create new df to predict on new values of x
new_dat <- data.frame(x = seq(from = 0, to = 1200, length.out = 100))
# Predict new fitted values and SE's on logodds scale
pred_logodds <- predict(mod_logistic, newdata = new_dat, type = "link", se = TRUE)
new_dat <- cbind(new_dat, pred_logodds)
# Create new df of predictions
predictions <- new_dat |> 
  rename(pred_logodds_est = fit) |> 
  mutate(pred_logodds_LL = pred_logodds_est - (1.96 * se.fit),
         pred_logodds_UL = pred_logodds_est + (1.96 * se.fit)) |> 
  select(-c(se.fit, residual.scale))
# Predict new fitted values and SE's on odds scale
predictions <- predictions |> 
  mutate(pred_odds_est = exp(pred_logodds_est),
         pred_odds_LL = exp(pred_logodds_LL),
         pred_odds_UL = exp(pred_logodds_UL))
# Predict new fitted values and SE's on probability scale
pred_probs <- predict(mod_logistic, newdata = new_dat, type = "response", se = TRUE)
new_dat <- cbind(new_dat[1], pred_probs)
new_dat <- new_dat |> 
  mutate(pred_probs_LL = fit - (1.96 * se.fit),
         pred_probs_UL = fit + (1.96 * se.fit))
# Add predicted probs and CIs to predictions df
predictions <- cbind(predictions, 
                     pred_probs_est = new_dat$fit, 
                     pred_probs_LL = new_dat$pred_probs_LL,
                     pred_probs_UL = new_dat$pred_probs_UL)

# Plots ----
p1 <- ggplot(predictions, aes(x = x, y = pred_logodds_est)) + 
  geom_ribbon(aes(ymin = pred_logodds_LL, ymax = pred_logodds_UL), alpha = 0.2) + 
  geom_line(color = "cornflowerblue", linewidth = 1) +
  geom_point(data = df, aes(x = x, y = y), size = 2, alpha = 0.1) +
  annotate("text", x = 1150, y = 6, label = "log-odds", size = 10, color = "#E7B800") +
  scale_x_continuous(limits = c(0, 1200), breaks = seq(0, 1200, by = 100)) +
  scale_y_continuous(limits = c(-100, 100), breaks = seq(-100, 100, by = 2)) +
  coord_cartesian(xlim = c(0, 1200), ylim = c(-8, 8)) +
  geom_vline(xintercept = 600, color = "red", linetype = "dotted", linewidth = 0.6) +
  geom_vline(xintercept = 700, color = "red", linetype = "dotted", linewidth = 0.6) +
  geom_label(aes(x = 630, y = 2), label = "1.323", color = "red", size = 3) +
  geom_segment(aes(x = 615, y = 1.72, xend = 600, yend = 1.323, color = "red"), 
               linewidth = 0.2, alpha = 0.05, arrow = arrow(length = unit(0.01, "npc"))) +
  geom_label(aes(x = 670, y = -1), label = "-0.181", color = "red", size = 3) +
  geom_segment(aes(x = 686, y = -0.7, xend = 700, yend = -0.181, color = "red"), 
               linewidth = 0.2, alpha = 0.05, arrow = arrow(length = unit(0.01, "npc"))) +
  guides(color = "none") + 
  theme_bw(base_size = 20) +
  ylab("") +
  theme(axis.title.x = element_blank(), axis.text.x = element_blank())
p2 <- ggplot(predictions, aes(x = x, y = pred_odds_est)) + 
  geom_ribbon(aes(ymin = pred_odds_LL, ymax = pred_odds_UL), alpha = 0.2) + 
  geom_line(color = "cornflowerblue", linewidth = 1) +
  geom_point(data = df, aes(x = x, y = y), size = 2, alpha = 0.1) +
  annotate("text", x = 1150, y = 6000, label = "odds", size = 10, color = "#E7B800") +
  scale_x_continuous(limits = c(0, 1200), breaks = seq(0, 1200, by = 100)) +
  scale_y_continuous(limits = c(-20, 1000000), breaks = c(0, 1 ,2, 3, 4, 5, seq(1000, 1000000, by = 1000))) +
  coord_cartesian(xlim = c(0, 1200), ylim = c(0, 7000)) +
  geom_vline(xintercept = 600, color = "red", linetype = "dotted", linewidth = 0.6) +
  geom_vline(xintercept = 700, color = "red", linetype = "dotted", linewidth = 0.6) +
  geom_label(aes(x = 630, y = 4), label = "3.755", color = "red", size = 3) +
  geom_segment(aes(x = 615, y = 3.85, xend = 600, yend = 3.755, color = "red"), 
               linewidth = 0.2, alpha = 0.05, arrow = arrow(length = unit(0.01, "npc"))) +
  geom_label(aes(x = 730, y = 1.5), label = "0.834", color = "red", size = 3) +
  geom_segment(aes(x = 715, y = 1.34, xend = 700, yend = 0.834, color = "red"), 
               linewidth = 0.2, alpha = 0.05, arrow = arrow(length = unit(0.01, "npc"))) +
  guides(color = "none") + 
  ylab("") +
  theme_bw(base_size = 20) +
  theme(axis.title.x = element_blank(), axis.text.x = element_blank()) +
  geom_magnify(from = c(xmin = 500, xmax = 1000, ymin = 0, ymax = 5), 
               to = c(xmin = 480, xmax = 1000, ymin = 1000, ymax = 5000), shadow = T, axes = "xy")
p3 <- ggplot(predictions, aes(x = x, y = pred_probs_est)) + 
  geom_ribbon(aes(ymin = pred_probs_LL, ymax = pred_probs_UL), alpha = 0.2) + 
  geom_line(color = "cornflowerblue", linewidth = 1) +
  geom_point(data = df, aes(x = x, y = y), size = 2, alpha = 0.1) +
  annotate("text", x = 1150, y = 0.8, label = "probability", size = 10, color = "#E7B800") +
  scale_x_continuous(limits = c(0, 1200), breaks = seq(0, 1200, by = 100)) +
  scale_y_continuous(limits = c(0, 1), breaks = seq(0, 1, by = 0.2)) +
  geom_vline(xintercept = 600, color = "red", linetype = "dotted", linewidth = 0.6) +
  geom_vline(xintercept = 700, color = "red", linetype = "dotted", linewidth = 0.6) +
  ylab("") + xlab("Macular hole thickness") +
  theme_bw(base_size = 20)
ggarrange(p1, p2, p3, align = "v", ncol = 1)
ggsave("combined.pdf", width = 20, height = 16)
[combined.pdf](https://github.com/hughjonesd/ggmagnify/files/13574906/combined.pdf)

# preliminary data transformation
[IND_COLL_data_Rus.csv](https://github.com/hughjonesd/ggmagnify/files/14462368/IND_COLL_data_Rus.csv)

# load data for plotting
df_classes <- 
![Figure 1_example](https://github.com/hughjonesd/ggmagnify/assets/93224567/5bd7efcd-a04d-4550-a19c-7d844771440e)
read.csv("IND_COLL_data_Rus.csv")

# obtain geoSpatial data for Russia
gadm <- raster::getData('GADM', country='RUS', level=1)
gadm1 <- raster::getData('GADM', country='UKR', level=1)
gadm1 <- gadm1[gadm1$NAME_1=="Crimea", ]
gadm <- rbind(gadm, gadm1)

for(i in 1:length(gadm@polygons)){
  for(j in 1:length(gadm@polygons[[i]]@Polygons)){
    gadm@polygons[[i]]@Polygons[[j]]@coords[,1]<- sapply(gadm@polygons[[i]]@Polygons[[j]]@coords[,1], function(x){
      if(x < 0){
        x<-359.999+x
      }
      else{x}
    })
  }
}


data_map <- tidy(gadm)
data_map$id <- as.numeric(data_map$id)
data_map1 <- left_join(data_map, df_classes, by =c("id" = "n_map"))
data_map1$na <- as.factor(ifelse(is.na(data_map1$IND_COLL), 1, 0))
table(data_map1$na)
tapply(data_map1$id, data_map1$oblast, mean) 
colnames(data_map1 )

 Define the location of a magnified map of St. Petersburg
from_SPb <- c(xmin = 28, xmax = 32, ymin = 58, ymax = 62)
# Names xmin, xmax, ymin, ymax are optional:
to_SPb <- c(20, 40, 72, 80)

# Define the location of a magnified map of Moscow
from_M <- c(xmin = 36, xmax = 39, ymin = 55, ymax = 58)
to_M <- c(70, 85, 41, 47)

# create the cultural map of Russia
map_base <- ggplot(data_map1, aes(x = long, y = lat, group = group, colour = IND_COLL)) + 
  geom_polygon(data = data_map1, colour='NA', linewidth = 0.01,
               aes(x=long, y=lat, fill=IND_COLL)) +
  #geom_polygon_pattern(data = data_map1, colour='black', size = 0.1,
   #                    aes(x=long, y=lat, fill=IND_COLL, pattern = na), 
    #                   pattern_color = "grey20",
    #                   pattern_fill = "grey20",
     #                  #pattern_density = 0.1,
     #                  pattern_linetype = 2,
      #                 pattern_spacing = 0.01,
      #                 pattern_size = 0.005) +
  #scale_pattern_manual(values = c('none', 'stripe'), guide = "none") +
  scale_fill_gradient(low = "gray1", high = "gray95", na.value = "white")+
  geom_polygon(data = subset(data_map1, id == 38), colour='white', fill = 'white')+
  theme(axis.line=element_blank(),axis.text.x=element_blank(),
        axis.text.y=element_blank(),axis.ticks=element_blank(),
        axis.title.x=element_blank(),
        axis.title.y=element_blank(),
        legend.position = c(0.85, 0.1),
        legend.direction = "horizontal",
        legend.key.size = unit(0.3, 'cm'),
        legend.text = element_text(colour = "black", #brewer.pal(11, 'RdYlGn')[4:9], 
                                      size = 4),
        legend.box.margin = margin(t = 0, b = 0),
        plot.margin = margin(t = 0, 
                             r = 0,  
                             b = 0,  
                             l = 0),
        panel.grid = element_blank(),
        panel.spacing = unit(0, "cm", data = NULL),
        panel.background = element_blank())+
  labs(x=NULL, y = NULL, fill = NULL)

map_base <- map_base+
  geom_magnify(from = from_SPb, to = to_SPb, linewidth = 0.15) +
  geom_magnify(from = from_M, to = to_M, linewidth = 0.15) 

  

tiff('Figure 1.tiff', units="cm", width=13.5, height=8, res=2000)
map_base # +
  #geom_magnify(from = from_M, to = to_M, linewidth = 0.15)  

# x coordinates for legend end labels have to be tuned manually which is extremely annoying
grid::grid.draw(grid::textGrob("К", x = 0.788, y = 0.13, 
                               gp = gpar(fontsize = 8, fontface = "bold")))
grid::grid.draw(grid::textGrob("И", x = 0.9315, y = 0.13, 
                               gp = gpar(fontsize = 8, fontface = "bold")))
dev.off()

ggp <- ggplot(iris, aes(Sepal.Width, Sepal.Length, color = Species)) + geom_point()
ggp


# Attempt 1: inset_theme() + plot()
axes <- "bla"  # without this line - Error: object 'axes' not found
Inset = ggp + inset_theme(
  blank = inset_blanks(axes = axes),
  axes,
  margin = if (axes == "") 50 else 0  # swapping '0' and '50' around shows changes in 'Inset'
)
Inset

ggm <- ggp + 
  geom_magnify(
    from = c(3, 3.5, 4.5, 5),
    to = c(3.5, 4.25, 6, 7),
    axes = "xy",
    plot = Inset)
ggm  # main plot and inset margins are not affected by if statement in inset_theme()


# Attempt 2: inset_theme() alone
ggm <- ggp + 
  geom_magnify(
    from = c(3, 3.5, 4.5, 5),
    to = c(3.5, 4.25, 6, 7),
    axes = "xy") +
  inset_theme(
    blank = inset_blanks(axes = axes),
    axes,
    margin = if (axes == "") 50 else 0
  )
ggm  # only main plot's margins are adjusted, inset plot margins remain unchanged


# Attempt 3: as Attempt 1, but just theme()
Inset = ggp + theme(
  plot.margin = unit(c(0, 0, 0, 0), "mm")) # changing the values from '0' changes the margins as expected
Inset

ggm <- ggp + 
  geom_magnify(
    from = c(3, 3.5, 4.5, 5),
    to = c(3.5, 4.25, 6, 7), axes = "xy", plot = Inset)
ggm # neither the main plot's nor the inset's margins are adjusted


# Attempt 4: as Attempt 2,  but just theme()
ggm <- ggp + 
  geom_magnify(
    from = c(3, 3.5, 4.5, 5),
    to = c(3.5, 4.25, 6, 7), axes = "xy") +
  theme(
    plot.margin = unit(c(20, 20, 20, 20), "mm"))
ggm # only the main plot's margins are adjusted

  gpanel <- ggpubr::as_ggplot(cowplot::get_legend(ggplot_object))
  gzoom <- ggmagnify(gpanel ,xlim = c(0.8, 1), ylim = c(0.8, 1),
                     inset_xlim = c(0, 0.5), inset_ylim = c(0, 0.5), zoom=2,
                     shadow = TRUE)

cost_table_ncbi <- readxl::read_excel('data/Sequencing_Cost_Data_Table_May2022.xls')

cost_table_ncbi$Date <- as.Date(cost_table_ncbi$Date)

# Define a custom formatter function
dollar_format <- function(x) {
  #ifelse(x >= 1e6, paste0("$", round(x / 1e6, 2), "mio"), paste0("$", comma(x)))
  paste0("$", comma(x))
}

ggplot(data = cost_table_ncbi, aes(x = Date)) +
  geom_line(aes(y = `Cost per Genome`), color = "blue") +
  coord_cartesian(clip = "off") + 
  ggtitle("Cost for sequencing a single human whole genome") +
  theme(axis.title.y = element_blank(),
        plot.title = element_text(face = "bold"))+
  scale_y_continuous(labels = dollar_format) +
  geom_magnify(from=list(ymd("2012-01-31"), ymd("2022-10-31"), 0 , 1e4),
               to=list(ymd("2008-05-31"), ymd("2023-5-31"),2.5e7, 8.5e7 ),
               shadow = TRUE, axes = "xy", recompute = TRUE, expand =0,
               after_scale = list(scale_y_continuous = dollar_format))

ggp <- ggplot(dv, aes(Position, NegLogP, color = cut)) + 
  geom_point(color = "darkblue", alpha = 0.8, size = 0.8) +
  labs(title = "GWAS p-values for cognitive function",
       subtitle = "Davies et al. (2018).", y = "-log(p)")

ggpi <- ggplot(iris, aes(Sepal.Width, Sepal.Length, colour = Species)) +
    geom_point() + xlim(c(1.5, 6))

p <- ggpi + geom_magnify(aes(from = Species == "setosa" & Sepal.Length < 5), 
                    to = c(4, 6, 6, 7.5))

ggsave("filename.eps", p, device = cairo_pdf)