---
title: "Case Study: Diabetic Nephropathy Meta-Analysis"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Case Study: Diabetic Nephropathy Meta-Analysis}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

This vignette demonstrates a comprehensive meta-analysis workflow for
diabetic nephropathy (DN) gene expression datasets. We show how to
systematically search, filter, annotate, and visualize relevant GEO datasets
to prepare for downstream meta-analysis.

## Introduction

Diabetic nephropathy is a major complication of diabetes and a leading cause
of chronic kidney disease. Conducting a meta-analysis of gene expression
studies can help identify robust molecular signatures. This workflow
demonstrates how to use `geokit` to:

- Search GEO using multiple query strategies
- Build a customized metadata database
- Filter datasets by quality criteria
- Visualize dataset characteristics
- Prepare a curated list for downstream analysis

```{r setup}
library(geokit)
library(dplyr)
library(stringr)
library(ggplot2)
```

## 1. Multi-Strategy Dataset Discovery

Different researchers use varying terminology, so we query for both
"diabetic nephropathy" and "diabetic kidney disease" and merge the results.

```{r search_dn, cache = TRUE}
dn_search_terms <- c(
  "diabetic nephropathy[ALL] AND Homo sapiens[ORGN] AND GSE[ETYP]",
  "diabetic kidney disease[ALL] AND Homo sapiens[ORGN] AND GSE[ETYP]"
)

dn_gse_list <- lapply(dn_search_terms, geo_search)
dn_gse <- unique(dplyr::bind_rows(dn_gse_list))

cat(sprintf("Found %d unique GSE datasets\n", nrow(dn_gse)))
head(dn_gse[, 1:3])
```

## 2. Extract Sample Information

Extract the number of samples from the `Contains` field to enable filtering
by sample size. This count covers every sample in a series regardless of
group — cases, controls, and technical replicates — so treat it as an upper
bound; group sizes and balance must be verified from sample annotations
before a dataset is included in a meta-analysis.

```{r extract_samples}
dn_gse <- dn_gse |>
  dplyr::mutate(
    number_of_samples = str_match(
      Contains, "(\\d+) Samples?"
    )[, 2L, drop = TRUE],
    number_of_samples = as.integer(number_of_samples)
  )

summary(dn_gse$number_of_samples)
```

## 3. Filter by Quality Criteria

Apply filters for sample size and study type to focus on candidate
expression profiling studies. Note that this filter admits both microarray
and RNA-seq records ("expression profiling by high throughput sequencing");
intensity- and count-based studies cannot be pooled directly, so process
them separately or tighten the filter to arrays for a single integrated
analysis.

```{r filter_datasets}
dn_gse_filtered <- dn_gse |>
  dplyr::filter(
    # At least 6 samples for meaningful analysis
    number_of_samples >= 6,
    # Focus on expression profiling studies
    str_detect(Type, regex("expression profiling", ignore_case = TRUE)),
    # Exclude methylation studies (defensive: no current GEO Type string
    # passes the filters above, but this guards against future variants)
    !str_detect(Type, regex("methylation", ignore_case = TRUE)),
    # Exclude compound SuperSeries so each study has a single assay type
    !str_detect(Type, ";")
  )

cat(sprintf(
  "After filtering: %d datasets (%.1f%% of original)\n",
  nrow(dn_gse_filtered),
  100 * nrow(dn_gse_filtered) / nrow(dn_gse)
))
```

## 4. Build a Metadata Database

Fetch detailed metadata for the filtered datasets using `geo_meta()`. This
creates a local database for offline analysis and detailed inspection. The
download may take several minutes depending on the number of datasets; if you
have an NCBI API key, set it as the `ENTREZ_KEY` environment variable to
speed up searching.

```{r build_metadb, eval = FALSE}
dn_metadb_dir <- "dn_metadb"
dir.create(dn_metadb_dir, showWarnings = FALSE, recursive = TRUE)

dn_metadb <- geo_meta(
  dn_gse_filtered[["Series Accession"]],
  odir = dn_metadb_dir
)

# Save for future use
saveRDS(dn_metadb, file.path(dn_metadb_dir, "dn_metadb.rds"))
```

In a later session, reload the saved database and calculate the actual sample
counts from the sample identifier list:

```{r load_metadb, eval = FALSE}
dn_metadb <- readRDS(file.path(dn_metadb_dir, "dn_metadb.rds"))

dn_metadb <- dn_metadb |>
  dplyr::mutate(
    number_of_samples = lengths(
      strsplit(Series_sample_id, "; ", fixed = TRUE)
    )
  )
```

## 5. Annotate and Categorize Datasets

Categorize datasets by criteria relevant to diabetic nephropathy research,
searching the title, summary, and overall design fields of each series.
These keyword assignments are screening heuristics only — inspect each
candidate's design and sample annotations (tissue, cell type, group
definitions) before final inclusion, as keyword matches can also pick up
cell-line or animal studies.

```{r categorize, eval = FALSE}
dn_metadb <- dn_metadb |>
  dplyr::mutate(
    search_text = paste(Series_title, Series_summary, Series_overall_design),
    has_nephropathy = str_detect(
      search_text, regex("nephropath", ignore_case = TRUE)
    ),
    has_kidney = str_detect(
      search_text, regex("kidney|renal", ignore_case = TRUE)
    ),
    tissue_type = dplyr::case_when(
      str_detect(search_text, regex("glomerul", ignore_case = TRUE)) ~
        "glomerular",
      str_detect(search_text, regex("tubul", ignore_case = TRUE)) ~ "tubular",
      str_detect(search_text, regex("kidney|renal", ignore_case = TRUE)) ~
        "kidney",
      TRUE ~ "other"
    )
  )

table(dn_metadb$tissue_type)
```

## 6. Visualization: Dataset Characteristics

### Timeline of Dataset Submissions

GSE accession numbers are assigned chronologically, so they can serve as a
time proxy for when studies were submitted. Depending on the records, the
SOFT metadata may also provide explicit dates (e.g., `Series_submission_date`)
that you can parse instead.

```{r timeline, eval = FALSE, fig.width = 8, fig.height = 4}
dn_metadb <- dn_metadb |>
  dplyr::mutate(
    gse_number = as.integer(str_extract(Series_geo_accession, "\\d+"))
  )

ggplot(dn_metadb, aes(x = gse_number)) +
  geom_histogram(bins = 30, fill = "steelblue", color = "white") +
  labs(
    title = "Distribution of DN Datasets by GEO Accession Number",
    subtitle = "GSE numbers are assigned chronologically (lower = older)",
    x = "GSE Number",
    y = "Number of Datasets"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    plot.subtitle = element_text(hjust = 0.5, size = 9),
    panel.grid.minor = element_blank()
  )
```

### Sample Size Distribution

```{r sample_distribution, eval = FALSE, fig.width = 8, fig.height = 4}
ggplot(dn_metadb, aes(x = number_of_samples)) +
  geom_histogram(binwidth = 5, fill = "coral", color = "white") +
  geom_vline(
    xintercept = median(dn_metadb$number_of_samples, na.rm = TRUE),
    linetype = "dashed", color = "red", linewidth = 1
  ) +
  labs(
    title = "Distribution of Sample Sizes in DN Datasets",
    subtitle = sprintf(
      "Median: %s samples",
      format(median(dn_metadb$number_of_samples, na.rm = TRUE))
    ),
    x = "Number of Samples",
    y = "Number of Datasets"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    plot.subtitle = element_text(hjust = 0.5),
    panel.grid.minor = element_blank()
  )
```

### Platform Usage Trends

```{r platform_trends, eval = FALSE, fig.width = 10, fig.height = 5}
platform_summary <- dn_metadb |>
  dplyr::count(Series_platform_id, sort = TRUE) |>
  dplyr::slice_max(n, n = 10)

ggplot(platform_summary, aes(x = reorder(Series_platform_id, n), y = n)) +
  geom_col(fill = "darkgreen", alpha = 0.7) +
  geom_text(aes(label = n), hjust = -0.2, size = 3.5) +
  coord_flip() +
  labs(
    title = "Top 10 Platforms Used in DN Studies",
    x = "Platform ID",
    y = "Number of Studies"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    panel.grid.major.y = element_blank()
  )
```

## 7. Prepare a Curated Dataset List

Create a final curated list of datasets suitable for meta-analysis.

```{r curate_list, eval = FALSE}
dn_curated <- dn_metadb |>
  dplyr::filter(
    number_of_samples >= 10,
    has_nephropathy | has_kidney
  ) |>
  dplyr::arrange(desc(number_of_samples)) |>
  dplyr::select(
    Series_geo_accession,
    Series_title,
    number_of_samples,
    Series_platform_id,
    tissue_type
  )

head(dn_curated, 10)

write.csv(
  dn_curated,
  file.path(dn_metadb_dir, "dn_curated_datasets.csv"),
  row.names = FALSE
)
```

## 8. Next Steps for Meta-Analysis

With this curated list, you can proceed to:

1. **Download expression data**: Use `geo_matrix()` to download series matrix
   files
2. **Quality control**: Examine sample annotations and expression
   distributions
3. **Data integration**: Normalize and batch-correct across studies
4. **Differential expression**: Identify genes consistently dysregulated in
   DN
5. **Pathway analysis**: Investigate enriched biological processes

```{r download_example, eval = FALSE}
top_gse <- dn_curated$Series_geo_accession[1]
eset <- geo_matrix(top_gse, odir = tempdir())

print(eset)
dim(Biobase::exprs(eset))
```

Instead of manually inspecting the curated list, you can also hand it to an
LLM for interactive exploration with `geo_chat()`; see
`vignette("geo-chat")` for details.

```{r chat_example, eval = FALSE}
chat <- geo_chat(
  NULL,
  dn_curated,
  data_description = paste(
    "Curated diabetic nephropathy GEO Series selected for meta-analysis.",
    "number_of_samples is counted from the SOFT sample list."
  )
)
chat$app()
```

## Summary

This workflow demonstrated how to:

- Use multiple search strategies to comprehensively identify relevant
  datasets
- Build a customized metadata database with `geo_meta()`
- Apply quality filters for sample size and study type
- Visualize dataset characteristics over time
- Prepare a curated list for downstream meta-analysis

The systematic approach supports reproducibility and helps identify the most
suitable datasets for gene expression meta-analysis in diabetic nephropathy
research.

## Session Information
```{r}
sessionInfo()
```