#your_email <- "czj9zj@virginia.edu"  # Replace with your email for better API performance

# Fetch all works
#works_5111646381 <- get_author_works("A5111646381", email = your_email)  # 43 works
#works_5112039075 <- get_author_works("A5112039075", email = your_email)  # 13 works
#works_5028976637 <- get_author_works("A5028976637", email = your_email)  # 2 works
#works_5102177495 <- get_author_works("A5102177495", email = your_email)  # 1 work

#data_path <- "/standard/nsdpi_storage/people/czj9zj/extra_data"
# Save each author's works
#saveRDS(works_5111646381, file = file.path(data_path, "works_5111646381.rds"))
#saveRDS(works_5112039075, file = file.path(data_path, "works_5112039075.rds"))
#saveRDS(works_5028976637, file = file.path(data_path, "works_5028976637.rds"))
#saveRDS(works_5102177495, file = file.path(data_path, "works_5102177495.rds"))

# Assuming your nuclear_scientists dataset is already loaded as 'nuclear_scientists'

# Merge the datasets
#enhanced_nuclear_scientists <- nuclear_scientists %>%left_join(author_info, by = "Victim")

# Display the enhanced dataset structure
#cat("Original columns:", ncol(nuclear_scientists), "\n")
#cat("Enhanced columns:", ncol(enhanced_nuclear_scientists), "\n")
#cat("New columns added:", ncol(enhanced_nuclear_scientists) - ncol(nuclear_scientists), "\n")

# Show column names
#cat("\nNew columns added:\n")
#new_cols <- setdiff(names(enhanced_nuclear_scientists), names(nuclear_scientists))
#cat(paste("-", new_cols, collapse = "\n"))

#saveRDS(enhanced_nuclear_scientists, "/standard/nsdpi_storage/people/czj9zj/extra_data/nuclear_scientists.rds")

# Display only the rows with author information

#data_path <- "extra_data"
#nuclear_scientists <- readRDS(file.path(data_path, "nuclear_scientists.rds"))


#kable(enhanced_nuclear_scientists, caption = "Nuclear Scientists with OpenAlex Information")

# Step 1: Find common columns
#common_cols <- intersect(names(new_abbasi_works), names(abbasi_works))

# Step 2: Subset both data frames to only those columns
#abbasi_subset <- abbasi_works[, common_cols]
#new_abbasi_subset <- new_abbasi_works[, common_cols]

# Step 3: Combine
#new_abbasi_works <- rbind(new_abbasi_subset, abbasi_subset)

#majid_id <- "A5112039075"

#abbasi_id <- "A5103401909"
#masoud_id <- "A5111436477"
#new_abbasi_id <- "a5037724110"

#your_email <- "czj9zj@virginia.edu"  # Replace with your email for better API performance

# Fetch all works
#cat("Fetching works for author ID:", author_id, "\n")
#majid_works <- get_author_works(majid_id, email = your_email)
#abbasi_works <- get_author_works(abbasi_id, email = your_email)


#new_abbasi_works <- get_author_works(new_abbasi_id, email = your_email)

#masoud_works <- get_author_works(masoud_id, email = your_email)

#data_path <- "/standard/nsdpi_storage/people/czj9zj/extra_data"
# Save Majid Shahriari's works
#saveRDS(majid_works, file = file.path(data_path, "majid_works.rds"))

# Save Fereydoon Abbasi's works
#saveRDS(abbasi_works, file = file.path(data_path, "abbasi_works.rds"))

#saveRDS(new_abbasi_works, file = file.path(data_path, "new_abbasi_works.rds"))

# Save Masoud Ali-Mohammadi's works
#masoud_works <- masoud_works %>%mutate(id = sub(".*/", "", id))
#saveRDS(masoud_works, file = file.path(data_path, "masoud_works.rds"))

#new_abbasi_works <- new_abbasi_works[, 1:50]

# Remove 'https://openalex.org/' from the 'id' column
#new_abbasi_works$id <- gsub("https://openalex.org/", "", new_abbasi_works$id)

#abbasi_works <- abbasi_works[, 1:50]

# Remove 'https://openalex.org/' from the 'id' column

#abbasi_works$id <- gsub("https://openalex.org/", "", abbasi_works$id)

# View the cleaned dataset
#head(abbasi_works)

# Step 1: Initialize an empty vector to store cleaned topic IDs
#all_topic_ids <- c()

# Step 2: Loop through all items in new_abbasi_works$topics
#for (i in seq_along(new_abbasi_works$topics)) {
  # Get the current topic dataframe
#  topic_df <- new_abbasi_works$topics[[i]]
  
  # Check if it's a data frame and contains the 'id' column
#  if (is.data.frame(topic_df) && "id" %in% names(topic_df)) {
    # Remove 'https://openalex.org/' from IDs
#    cleaned_ids <- gsub("https://openalex.org/", "", topic_df$display_name)
    
    # Append to the list of all topic IDs
#    all_topic_ids <- c(all_topic_ids, cleaned_ids)
#  }
#}

# Step 3: Optional - remove duplicates
#topics_name <- unique(all_topic_ids)
#new_abbasi_keywords <- topic_summary[topic_summary$Topic %in% topics_name, ]
#saveRDS(new_abbasi_keywords, file = file.path(data_path, "new_abbasi_keywords.rds"))

data_path <- "/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic"
majid_works <- readRDS(file.path(data_path, "MShahriari_works.rds"))
abbasi_works <- readRDS(file.path(data_path, "FAbbasiDavani_works.rds"))
shahid_beheshti_university_cleanedLONG <-readRDS("/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data/shahid_beheshti_university_cleanedLONG.rds")

                                                  masoud_works <- readRDS(file.path("/standard/nsdpi_storage/people/czj9zj/extra_data", "masoud_works.rds"))

# Define your path if needed
data_path <- "/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data"

# Read the keyword files (adjust file names if needed)
#majid_keywords <- readRDS(file.path(data_path, "majid_keywords.rds"))
#abbasi_keywords <- readRDS(file.path(data_path, "abbasi_keywords.rds"))

majid_keywords <- readRDS(file.path(data_path, "majid_keywords_v1.rds"))
new_abbasi_keywords <- readRDS(file.path(data_path, "abbasi_keywords_v1.rds"))

# Step 1: Extract unique topics from the 'Topic' column
#abbasi_topics <- unique(new_abbasi_keywords$Topic)
#majid_topics <- unique(new_majid_keywords$Topic)

# Step 2: Identify topic categories
#shared_topics <- intersect(abbasi_topics, majid_topics)
#abbasi_only_topics <- setdiff(abbasi_topics, shared_topics)
#majid_only_topics <- setdiff(majid_topics, shared_topics)

# Step 3: Padding helper
#pad_column <- function(x, n) {
#  length(x) <- n
#  x[is.na(x)] <- "---"
#  return(x)
#}

# Step 4: Create the comparison table
#max_len <- max(length(majid_only_topics), length(shared_topics), length(abbasi_only_topics))

#topic_comparison <- tibble(
#  "Majid Shahriari Only" = pad_column(majid_only_topics, max_len),
#  "Shared Topics" = pad_column(shared_topics, max_len),
#  "Fereydoon Abbasi Only" = pad_column(abbasi_only_topics, max_len)
#)

# Step 5: Generate the styled table
#table_html <- kbl(
#  topic_comparison,
#  caption = "Topic Comparison: Majid Shahriari vs. Fereydoon Abbasi"
#) %>%
#  kable_styling(full_width = FALSE)

# Step 6: Save as HTML
#dir.create("plots", showWarnings = FALSE)
#save_kable(table_html, file = "plots/topic_comparison.html")

# Step 7: Save as PNG
#webshot("plots/topic_comparison.html", "plots/topic_comparison.png", zoom = 2)

length(unique(works_393$id)) == nrow(works_393)

length(unique(works_403$id)) == nrow(works_403)

duplicated_ids <- works_403$id[duplicated(works_403$id)]

# Show which rows are duplicated
works_403 %>% 
  filter(id %in% duplicated_ids)

works_403 <- works_403 %>%
  filter(!duplicated(id) & !duplicated(id, fromLast = TRUE))

  domain_colors <- c(
    "Physical Sciences" = "#1e3a8a",    # Deep blue
    "Social Sciences"   = "#7c2d12",    # Deep brown/orange
    "Life Sciences"     = "#14532d",    # Deep green  
    "Health Sciences"   = "#7f1d1d"     # Deep red
  )
  

# Create the plot and assign it to domain_pct_plot
domain_pct_plot <- works_393 %>%
  filter(!is.na(top_domain)) %>%
  count(top_domain, sort = TRUE) %>%
  mutate(pct = n / sum(n)) %>%
  ggplot(aes(x = reorder(top_domain, pct), y = pct, fill = top_domain)) +
  geom_col() +
  geom_text(aes(label = percent(pct, accuracy = 0.1)),
            hjust = 1.3, color = "white", size = 4.2) +
  coord_flip() +
  scale_y_continuous(labels = percent_format(accuracy = 0.1)) +
  scale_fill_manual(values = domain_colors) +
  labs(
    title = "Research Domains at Shahid Beheshti University",
    x = "Domain",
    y = "% of Total Works",
    fill = "Domain"
  ) +
  theme_minimal() +
  theme(
    legend.position = "none",
    axis.text = element_text(size = 12),
    axis.title = element_text(size = 13, face = "bold"),
    plot.title = element_text(size = 16, face = "bold", hjust = 0.5)
  )

# Now save it
#ggsave("plots/domain_percent_bar.png", plot = domain_pct_plot, width = 10, height = 6, dpi = 320)

topics_data <- works_393 %>%
  filter(!is.na(top_domain), !is.na(top_field)) %>%
  group_by(domain.display_name = top_domain, field.display_name = top_field) %>%
  summarise(works_count = n(), .groups = "drop")

treemap_plot <- create_complete_names_treemap(topics_data)

#ggsave("plots/domain_treemap.png", treemap_plot, width = 12, height = 8, dpi = 320)

# Count unique fields
n_fields <- works_393 %>%
  filter(!is.na(top_field)) %>%
  summarise(n = n_distinct(top_field)) %>%
  pull(n)

# Get field-level counts and domains
field_data <- works_393 %>%
  filter(!is.na(top_field), !is.na(top_domain)) %>%
  count(top_field, top_domain, sort = TRUE)

# Plot all fields
ggplot(field_data, aes(x = reorder(top_field, n), y = n, fill = top_domain)) +
  geom_col() +
  coord_flip() +
  scale_fill_manual(values = domain_colors) +
  labs(
    title = "Distribution of Research Fields at Shahid Beheshti University",
    subtitle = paste("All", n_fields, "fields colored by domain"),
    x = "Field",
    y = "Number of Works",
    fill = "Domain"
  ) +
  theme_minimal() +
  theme(
    axis.text = element_text(size = 11),
    plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
    plot.subtitle = element_text(size = 11, hjust = 0.5)
  )

library(tidyverse)

library(forcats)

# Identify top 5 fields by publication count
top_5_fields <- works_393 %>%
  filter(!is.na(top_field), !is.na(publication_year)) %>%
  count(publication_year, top_field) %>%
  group_by(top_field) %>%
  mutate(total = sum(n)) %>%
  ungroup() %>%
  group_by(top_field) %>%
  slice_max(order_by = total, n = 1, with_ties = FALSE) %>%
  ungroup() %>%
  slice_max(order_by = total, n = 5) %>%
  pull(top_field)

# Filter and plot
plot <- works_393 %>%
  filter(top_field %in% top_5_fields, !is.na(publication_year)) %>%
  count(publication_year, top_field) %>%
  ggplot(aes(x = publication_year, y = n, fill = top_field)) +
  geom_col(show.legend = FALSE) +
  facet_wrap(~ top_field, scales = "free_y", ncol = 2) +
  scale_x_continuous(limits = c(2000, 2025), breaks = seq(2000, 2025, 2)) +
  theme_minimal() +
  labs(
    title = "Evolution of Top 5 Research Fields Over Time",
    subtitle = "Shahid Beheshti University – by publication count",
    x = "Publication Year",
    y = "Number of Works"
  )

# Save the plot
#ggsave("plots/top_fields_evolution.png", plot = plot, width = 10, height = 6, dpi = 300)

library(tidyverse)
library(forcats)

# Identify top 5 subfields by publication count
top_5_subfields <- works_393 %>%
  filter(!is.na(top_subfield), !is.na(publication_year)) %>%
  count(publication_year, top_subfield) %>%
  group_by(top_subfield) %>%
  mutate(total = sum(n)) %>%
  ungroup() %>%
  group_by(top_subfield) %>%
  slice_max(order_by = total, n = 1, with_ties = FALSE) %>%
  ungroup() %>%
  slice_max(order_by = total, n = 5) %>%
  pull(top_subfield)

# Filter and plot
plot <- works_393 %>%
  filter(top_subfield %in% top_5_subfields, !is.na(publication_year)) %>%
  count(publication_year, top_subfield) %>%
  ggplot(aes(x = publication_year, y = n, fill = top_subfield)) +
  geom_col(show.legend = FALSE) +
  facet_wrap(~ top_subfield, scales = "free_y", ncol = 2) +
  scale_x_continuous(limits = c(2000, 2025), breaks = seq(2000, 2025, 2)) +
  theme_minimal() +
  labs(
    title = "Evolution of Top 5 Research Subfields Over Time",
    subtitle = "Shahid Beheshti University – by publication count",
    x = "Publication Year",
    y = "Number of Works"
  )

plot

# Save the plot
#ggsave("plots/top_subfields_evolution.png", plot = plot, width = 10, height = 6, dpi = 300)

library(tidyverse)
library(broom)  # for tidy regression results

# Step 1: Filter and prepare data
field_growth_slopes <- works_393 %>%
  filter(!is.na(top_field), !is.na(publication_year)) %>%
  count(publication_year, top_field) %>%
  group_by(top_field) %>%
  filter(n() >= 5) %>%  # Ensure enough data points for a slope
  nest() %>%
  mutate(
    model = map(data, ~lm(n ~ publication_year, data = .x)),
    slope = map_dbl(model, ~coef(.x)[["publication_year"]])
  ) %>%
  ungroup() %>%
  arrange(desc(slope)) %>%
  slice_head(n = 5)

# Extract top 5 fastest-growing fields
top_5_growing_fields <- field_growth_slopes$top_field

# Filter and plot the fastest-growing fields
plot <- works_393 %>%
  filter(top_field %in% top_5_growing_fields, !is.na(publication_year)) %>%
  count(publication_year, top_field) %>%
  ggplot(aes(x = publication_year, y = n, fill = top_field)) +
  geom_col(show.legend = FALSE) +
  facet_wrap(~ top_field, scales = "free_y", ncol = 2) +
  scale_x_continuous(limits = c(2000, 2025), breaks = seq(2000, 2025, 2)) +
  theme_minimal() +
  labs(
    title = "Fastest-Growing Research Fields Over Time",
    subtitle = "Top 5 fields at Shahid Beheshti University based on publication growth rate",
    x = "Publication Year",
    y = "Number of Works"
  )

plot

#ggsave("plots/sbu_fastest_growing_fields.png", plot = plot, width = 10, height = 6, dpi = 300)

library(tidyverse)
library(scales)
library(lubridate)

# Ensure publication dates are parsed
works_393 <- works_393 %>%
  mutate(pub_month = floor_date(publication_date, "month"))

# Filter after year 2000 only
monthly_counts <- works_393 %>%
  filter(top_subfield %in% top_5_subfields, pub_month >= as.Date("2000-01-01")) %>%
  count(pub_month, top_subfield)

# Normalize to each subfield's first available month (after 2000)
monthly_growth <- monthly_counts %>%
  group_by(top_subfield) %>%
  arrange(pub_month) %>%
  mutate(
    base_count = first(n),
    growth_index = (n / base_count) * 100
  ) %>%
  ungroup()



dark_colors <- c(
  "Biomedical Engineering" = "#102542",   # Dark navy
  "Electrical and Electronic Engineering" = "#A84400",  # Deep indigo
  "Materials Chemistry" = "darkgreen",      # Dark red (true crimson)
  "Molecular Biology" = "darkred",        # Deep cyan-teal
  "Organic Chemistry" = "#4B0082"         # Dark violet/indigo
)


ggplot(monthly_growth, aes(x = pub_month, y = growth_index, color = top_subfield)) +
  geom_line(linewidth = 1, alpha = 0.8) +
  scale_color_manual(values = dark_colors) +
  scale_y_continuous(labels = percent_format(scale = 1)) +
  scale_x_date(date_breaks = "2 years", date_labels = "%Y", limits = c(as.Date("2000-01-01"), NA)) +
  labs(
    title = "Monthly Publication Growth in Top Subfields at Shahid Beheshti University",
    subtitle = "Indexed to each subfield's first recorded month after 2000 (Index = 100%)",
    x = "Publication Date",
    y = "Relative Publication Index (%)",
    color = "Subfield"
  ) +
  theme_minimal(base_size = 12)

# Save the plot
ggsave("plots/monthly_subfield_growth.png", width = 10, height = 6, dpi = 300)

#ggsave("plots/subfield_growth_slopechart.png", width = 10, height = 6, dpi = 300)

monthly_growth_sbu <- monthly_growth %>%
  mutate(university = "Shahid Beheshti University")

works_403 <- works_403 %>%
  mutate(pub_month = floor_date(publication_date, "month"))

# Filter after year 2000 only
monthly_counts_ut <- works_403 %>%
  filter(top_subfield %in% top_5_subfields, pub_month >= as.Date("2000-01-01")) %>%
  count(pub_month, top_subfield)

# Normalize to each subfield's first available month (after 2000)
monthly_growth_ut <- monthly_counts %>%
  group_by(top_subfield) %>%
  arrange(pub_month) %>%
  mutate(
    base_count = first(n),
    growth_index = (n / base_count) * 100
  ) %>%
  ungroup()


monthly_growth_ut <- monthly_growth_ut %>%
  mutate(university = "University of Tehran")
shared_colors <- c(
  "Biomedical Engineering" = "#102542",   # Dark navy
  "Electrical and Electronic Engineering" = "#A84400",  # Burnt orange
  "Molecular Biology" = "#8B0000",        # Dark red
  "Materials Chemistry" = "#4B0082"       # Dark violet
)

sbu_only_colors <- c(
  "Organic Chemistry" = "#3E3E6B"         # Deep steel blue
)

ut_only_colors <- c(
  "Mechanical Engineering" = "darkgreen"    # Dark teal
)

dark_colors_combined <- c(shared_colors, sbu_only_colors, ut_only_colors)

# Using new_abbasi_keywords instead of just abbasi_keywords
# Subset works in Fereydoon Abbasi’s areas
data_path <- "/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data"

abbasi_keywords_v1 <- readRDS(file.path(data_path, "abbasi_keywords_v1.rds"))


sbu_abbasi_works <- shahid_beheshti_university_cleanedLONG %>%
  filter(
    topic %in% abbasi_keywords_v1$Topic
  )

# Works NOT in Abbasi’s areas
sbu_non_abbasi_works <- shahid_beheshti_university_cleanedLONG %>%
  filter(
    !(topic %in% abbasi_keywords_v1$Topic)
  )


sbu_abbasi_works %>%
  head(15) %>%
  kbl(caption = "Subset of Works from Shahid Beheshti University Matching Fereydoon Abbasi’s Research Areas") %>%
  kable_styling(latex_options = c("hold_position", "scale_down"), font_size = 10)

library(dplyr)

data_path <- "/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data"

masoud_keywords_v1 <- readRDS(file.path(data_path, "masoud_keywords_v1.rds"))

base_path <- "/standard/nsdpi_storage/people/czj9zj/temp_progress_new"

# university_of_tehran_long <- readRDS(file.path(base_path, "works_batch_403.rds")) %>%
#   filter(type != "dataset", type != "erratum", type != "retraction") %>%
#   select(id, title, publication_date, topics, type) %>%
#   # keep rows where 'topics' exists and has at least one row
#   filter(purrr::map_lgl(topics, ~ !is.null(.x) && nrow(.x) > 0)) %>%
#   mutate(
#     topics = map(topics, ~ filter(.x, type == "topic")),
#     topics = map(topics, ~ select(.x, display_name, score))
#   ) %>%
#   unnest(cols = topics) %>%
#   rename(topic = display_name) %>%
#   mutate(
#     across(where(is.character), ~ gsub("https://openalex.org/", "", .x, fixed = TRUE)),
#     institution = "University of Tehran"
#   )


#saveRDS(university_of_tehran_long, file.path(base_path, "university_of_tehran_long.rds"))

university_of_tehran_long <- readRDS(file.path(base_path, "university_of_tehran_long.rds"))

# Works in Masoud's subfields
# Works in Masoud’s areas
sbu_masoud_works <- university_of_tehran_long %>%
  filter(
    topic %in% masoud_keywords_v1$Topic
  )

# Works NOT in Masoud’s areas
sbu_non_masoud_works <- university_of_tehran_long %>%
  filter(
    !(topic %in% masoud_keywords_v1$Topic)
  )

# Save both datasets as RDS
saveRDS(sbu_masoud_works, file = file.path("/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data", "sbu_masoud_works.rds"))
saveRDS(sbu_non_masoud_works, file = file.path("/standard/nsdpi_storage/people/czj9zj/DiD_data_3topic/extra_data", "sbu_non_masoud_works.rds"))