inst/tinytest/test_integration.R

In Rduckhts: 'DuckDB' High Throughput Sequencing File Formats Reader Extension

# Integration tests for Rduckhts package - shows RBCFTools-like table creation
library(tinytest)
library(DBI)

# Test table creation pattern (like RBCFTools)
test_table_creation <- function() {
  # Create an in-memory DuckDB connection with the bundled extension loaded.
  con <- rduckhts_connect()
  on.exit(dbDisconnect(con, shutdown = TRUE), add = TRUE)

  runtime_probe <- dbGetQuery(
    con,
    paste(
      "SELECT",
      "duckhts_duckdb_type_supported('VARCHAR') AS has_varchar,",
      "duckhts_duckdb_type_supported('variant') = duckhts_duckdb_supports_variant() AS variant_consistent,",
      "duckhts_duckdb_type_supported('geometry') = duckhts_duckdb_supports_geometry() AS geometry_consistent,",
      "NOT duckhts_duckdb_type_supported('__duckhts_no_such_type__') AS missing_is_false"
    )
  )
  expect_true(runtime_probe$has_varchar)
  expect_true(runtime_probe$variant_consistent)
  expect_true(runtime_probe$geometry_consistent)
  expect_true(runtime_probe$missing_is_false)

  bcf_path <- system.file("extdata", "vcf_file.bcf", package = "Rduckhts")
  bcf_index_path <- system.file(
    "extdata",
    "vcf_file.bcf.csi",
    package = "Rduckhts"
  )
  formatcols_vcf_path <- system.file(
    "extdata",
    "formatcols.vcf.gz",
    package = "Rduckhts"
  )
  formatcols_vcf_index_path <- system.file(
    "extdata",
    "formatcols.vcf.gz.csi",
    package = "Rduckhts"
  )
  bam_path <- system.file("extdata", "range.bam", package = "Rduckhts")
  bam_index_path <- system.file(
    "extdata",
    "range.bam.bai",
    package = "Rduckhts"
  )
  fasta_path <- system.file("extdata", "ce.fa", package = "Rduckhts")
  fastq_r1 <- system.file("extdata", "r1.fq", package = "Rduckhts")
  fastq_r2 <- system.file("extdata", "r2.fq", package = "Rduckhts")
  gff_path <- system.file("extdata", "gff_file.gff.gz", package = "Rduckhts")
  gff_index_path <- system.file(
    "extdata",
    "gff_file.gff.gz.tbi",
    package = "Rduckhts"
  )
  gff_strict_valid_path <- system.file(
    "extdata",
    "gff_strict_valid.gff3",
    package = "Rduckhts"
  )
  gff_strict_invalid_path <- system.file(
    "extdata",
    "gff_strict_invalid.gff3",
    package = "Rduckhts"
  )
  gff_strict_invalid_attr_path <- system.file(
    "extdata",
    "gff_strict_invalid_attr.gff3",
    package = "Rduckhts"
  )
  gff_strict_extra_field_path <- system.file(
    "extdata",
    "gff_strict_extra_field.gff3",
    package = "Rduckhts"
  )
  gff_strict_invalid_end_path <- system.file(
    "extdata",
    "gff_strict_invalid_end.gff3",
    package = "Rduckhts"
  )
  gff_attrs_path <- system.file(
    "extdata",
    "gff_attrs.gff3",
    package = "Rduckhts"
  )
  gtf_attrs_path <- system.file(
    "extdata",
    "gtf_attrs.gtf",
    package = "Rduckhts"
  )
  tabix_path <- system.file("extdata", "rg.sam.gz", package = "Rduckhts")
  header_tabix_path <- system.file(
    "extdata",
    "header_tabix.tsv.gz",
    package = "Rduckhts"
  )
  header_tabix_index_path <- system.file(
    "extdata",
    "header_tabix.tsv.gz.tbi",
    package = "Rduckhts"
  )
  meta_tabix_path <- system.file(
    "extdata",
    "meta_tabix.tsv.gz",
    package = "Rduckhts"
  )
  vep_path <- system.file("extdata", "test_vep.vcf", package = "Rduckhts")

  expect_true(file.exists(bcf_path))
  expect_true(file.exists(bcf_index_path))
  expect_true(file.exists(formatcols_vcf_path))
  expect_true(file.exists(formatcols_vcf_index_path))
  expect_true(file.exists(bam_path))
  expect_true(file.exists(bam_index_path))
  expect_true(file.exists(fasta_path))
  expect_true(file.exists(fastq_r1))
  expect_true(file.exists(fastq_r2))
  expect_true(file.exists(gff_path))
  expect_true(file.exists(gff_index_path))
  expect_true(file.exists(gff_strict_valid_path))
  expect_true(file.exists(gff_strict_invalid_path))
  expect_true(file.exists(gff_strict_invalid_attr_path))
  expect_true(file.exists(gff_strict_extra_field_path))
  expect_true(file.exists(gff_strict_invalid_end_path))
  expect_true(file.exists(gff_attrs_path))
  expect_true(file.exists(gtf_attrs_path))
  expect_true(file.exists(tabix_path))
  expect_true(file.exists(header_tabix_path))
  expect_true(file.exists(header_tabix_index_path))
  expect_true(file.exists(meta_tabix_path))
  expect_true(file.exists(vep_path))

  fasta_index_path <- tempfile("duckhts_fasta_", fileext = ".fai")
  expect_silent(rduckhts_fasta_index(
    con,
    fasta_path,
    index_path = fasta_index_path
  ))
  expect_true(file.exists(fasta_index_path))

  expect_silent(rduckhts_bcf(con, "variants", bcf_path, overwrite = TRUE))
  expect_silent(rduckhts_bcf(
    con,
    "variants_idx",
    bcf_path,
    region = "1:3000150-3000151",
    index_path = bcf_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_bcf(
    con,
    "variants_multi_idx",
    bcf_path,
    region = "1:3000150-3000151,1:3000150-3000151",
    index_path = bcf_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_bam(con, "reads", bam_path, overwrite = TRUE))
  expect_silent(rduckhts_bam(
    con,
    "reads_idx",
    bam_path,
    region = "CHROMOSOME_I:1-1000",
    index_path = bam_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_fasta(con, "sequences", fasta_path, overwrite = TRUE))
  expect_silent(rduckhts_fasta(
    con,
    "sequences_region",
    fasta_path,
    region = "CHROMOSOME_I:1-10",
    index_path = fasta_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_fastq(
    con,
    "fastq_reads",
    fastq_r1,
    mate_path = fastq_r2,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_gff(
    con,
    "annotations",
    gff_path,
    attributes_map = TRUE,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_gff(
    con,
    "annotations_idx",
    gff_path,
    region = "X:2934816-2935190",
    index_path = gff_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_gff(
    con,
    "annotations_multi_idx",
    gff_path,
    region = "X:2934816-2935190,X:2934816-2935190",
    index_path = gff_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_gff(
    con,
    "annotations_strict",
    gff_strict_valid_path,
    strict = TRUE,
    overwrite = TRUE
  ))
  expect_equal(
    DBI::dbGetQuery(
      con,
      "SELECT count(*) AS n, count(start) AS n_start, count(\"end\") AS n_end FROM annotations_strict"
    ),
    data.frame(n = 3, n_start = 2, n_end = 2)
  )
  expect_error(
    rduckhts_gff(
      con,
      "annotations_strict_bad",
      gff_strict_invalid_path,
      strict = TRUE,
      overwrite = TRUE
    ),
    pattern = "InvalidCoordinate"
  )
  expect_error(
    rduckhts_gff(
      con,
      "annotations_strict_bad_attr",
      gff_strict_invalid_attr_path,
      strict = TRUE,
      overwrite = TRUE
    ),
    pattern = "InvalidAttribute"
  )
  expect_error(
    rduckhts_gff(
      con,
      "annotations_strict_extra_field",
      gff_strict_extra_field_path,
      strict = TRUE,
      overwrite = TRUE
    ),
    pattern = "TooManyFields"
  )
  expect_error(
    rduckhts_gff(
      con,
      "annotations_strict_bad_end",
      gff_strict_invalid_end_path,
      strict = TRUE,
      overwrite = TRUE
    ),
    pattern = "InvalidCoordinate"
  )
  expect_silent(rduckhts_gff(
    con,
    "annotations_attrs",
    gff_attrs_path,
    attributes_list = TRUE,
    attributes_pairs = TRUE,
    overwrite = TRUE
  ))
  expect_equal(
    DBI::dbGetQuery(
      con,
      paste(
        "SELECT list_count(map_extract_value(attributes_list, 'Dbxref')) AS n_dbxref,",
        "list_extract(map_extract_value(attributes_list, 'Note'), 1) AS note,",
        "list_count(attributes_pairs) AS n_pairs FROM annotations_attrs"
      )
    ),
    data.frame(n_dbxref = 2, note = "hello world", n_pairs = 6)
  )
  expect_silent(rduckhts_gtf(
    con,
    "gtf_attrs",
    gtf_attrs_path,
    attributes_list = TRUE,
    attributes_pairs = TRUE,
    overwrite = TRUE
  ))
  expect_equal(
    DBI::dbGetQuery(
      con,
      paste(
        "SELECT list_extract(map_extract_value(attributes_list, 'note'), 1) AS note,",
        "list_count(attributes_pairs) AS n_pairs FROM gtf_attrs"
      )
    ),
    data.frame(note = "weird; semi", n_pairs = 3)
  )
  expect_silent(rduckhts_tabix(con, "tabix_data", tabix_path, overwrite = TRUE))
  expect_silent(rduckhts_tabix(
    con,
    "tabix_idx",
    gff_path,
    region = "X:2934816-2935190",
    index_path = gff_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_tabix(
    con,
    "tabix_multi_idx",
    gff_path,
    region = "X:2934816-2935190,X:2934816-2935190",
    index_path = gff_index_path,
    overwrite = TRUE
  ))
  expect_silent(rduckhts_bcf(con, "vep_variants", vep_path, overwrite = TRUE))
  expect_true(DBI::dbExistsTable(con, "variants_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM variants_idx")$n[1] == 2
  )
  expect_true(DBI::dbExistsTable(con, "variants_multi_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM variants_multi_idx")$n[1] ==
      2
  )
  expect_true(DBI::dbExistsTable(con, "reads_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM reads_idx")$n[1] == 2
  )
  expect_true(DBI::dbExistsTable(con, "tabix_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM tabix_idx")$n[1] == 4
  )
  expect_true(DBI::dbExistsTable(con, "tabix_multi_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM tabix_multi_idx")$n[1] == 4
  )
  expect_true(DBI::dbExistsTable(con, "annotations_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM annotations_idx")$n[1] == 4
  )
  expect_true(DBI::dbExistsTable(con, "annotations_multi_idx"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM annotations_multi_idx")$n[
      1
    ] ==
      4
  )
  expect_true(DBI::dbExistsTable(con, "sequences_region"))
  expect_true(
    DBI::dbGetQuery(con, "SELECT count(*) AS n FROM sequences_region")$n[1] == 1
  )
  expect_true(
    DBI::dbGetQuery(
      con,
      "SELECT length(SEQUENCE) AS n FROM sequences_region"
    )$n[1] ==
      10
  )

  bcf_count_only <- DBI::dbGetQuery(
    con,
    sprintf(
      paste(
        "SELECT COUNT(*) AS n FROM read_bcf(",
        "'%s', tidy_format := true, index_path := '%s')"
      ),
      bcf_path,
      bcf_index_path
    )
  )
  expect_equal(bcf_count_only$n[1], 30)

  gff_count_only <- DBI::dbGetQuery(
    con,
    sprintf(
      paste(
        "SELECT COUNT(*) AS n FROM read_gff(",
        "'%s', index_path := '%s')"
      ),
      gff_path,
      gff_index_path
    )
  )
  expect_equal(gff_count_only$n[1], 62)

  tabix_header_count_only <- DBI::dbGetQuery(
    con,
    sprintf(
      paste(
        "SELECT COUNT(*) AS n FROM read_tabix(",
        "'%s', header := TRUE, index_path := '%s')"
      ),
      header_tabix_path,
      header_tabix_index_path
    )
  )
  expect_equal(tabix_header_count_only$n[1], 2)

  header_meta <- rduckhts_hts_header(con, bcf_path)
  expect_true(nrow(header_meta) > 0)
  expect_true("record_type" %in% names(header_meta))
  header_raw <- rduckhts_hts_header(con, bcf_path, mode = "raw")
  expect_true(nrow(header_raw) > 0)
  expect_true("raw" %in% names(header_raw))

  index_meta <- rduckhts_hts_index(con, bcf_path, index_path = bcf_index_path)
  expect_true(nrow(index_meta) > 0)
  expect_true("index_type" %in% names(index_meta))
  index_spans <- rduckhts_hts_index_spans(
    con,
    formatcols_vcf_path,
    index_path = formatcols_vcf_index_path
  )
  expect_equal(nrow(index_spans), 1)
  expect_true("chunk_beg_vo" %in% names(index_spans))
  expect_equal(index_spans$seqname[[1]], "1")
  expect_equal(index_spans$tid[[1]], 0)
  expect_equal(index_spans$bin[[1]], 2396745)
  expect_equal(index_spans$chunk_beg_vo[[1]], 20381696)
  expect_equal(index_spans$chunk_end_vo[[1]], 23789568)
  expect_equal(index_spans$seq_start[[1]], 1)
  expect_equal(index_spans$seq_end[[1]], 16384)
  index_raw <- rduckhts_hts_index_raw(
    con,
    formatcols_vcf_path,
    index_path = formatcols_vcf_index_path
  )
  expect_true(nrow(index_raw) == 1)
  expect_true("raw" %in% names(index_raw))
  expect_true(DBI::dbExistsTable(con, "annotations"))
  if (DBI::dbExistsTable(con, "annotations")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM annotations LIMIT 1"))

    # Check MAP types specifically
    schema <- DBI::dbGetQuery(con, "DESCRIBE annotations")
    message(
      "GFF schema with attributes_map=TRUE:\n",
      paste(capture.output(schema), collapse = "\n")
    )

    # Check what type the attributes column has
    attr_type <- DBI::dbGetQuery(
      con,
      "SELECT typeof(attributes) as attr_type FROM annotations LIMIT 1"
    )$attr_type[1]
    message("Attributes column type: ", attr_type)

    # Get sample to see MAP structure
    sample <- DBI::dbGetQuery(con, "SELECT attributes FROM annotations LIMIT 1")
    message(
      "Sample MAP content:\n",
      paste(capture.output(sample), collapse = "\n")
    )
  }

  expect_true(DBI::dbExistsTable(con, "reads"))
  if (DBI::dbExistsTable(con, "reads")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM reads LIMIT 1"))
  }

  expect_true(DBI::dbExistsTable(con, "sequences"))
  if (DBI::dbExistsTable(con, "sequences")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM sequences LIMIT 1"))
  }

  expect_true(DBI::dbExistsTable(con, "fastq_reads"))
  if (DBI::dbExistsTable(con, "fastq_reads")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM fastq_reads LIMIT 1"))
  }

  expect_true(DBI::dbExistsTable(con, "annotations"))
  if (DBI::dbExistsTable(con, "annotations")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM annotations LIMIT 1"))
  }

  expect_true(DBI::dbExistsTable(con, "tabix_data"))
  if (DBI::dbExistsTable(con, "tabix_data")) {
    expect_silent(DBI::dbGetQuery(con, "SELECT * FROM tabix_data LIMIT 1"))
  }

  expect_silent(rduckhts_tabix(
    con,
    "tabix_header",
    header_tabix_path,
    header = TRUE,
    overwrite = TRUE
  ))
  tabix_header_cols <- DBI::dbGetQuery(con, "PRAGMA table_info('tabix_header')")
  expect_equal(tabix_header_cols$name, c("chrom", "pos", "value"))

  expect_silent(rduckhts_tabix(
    con,
    "tabix_named",
    meta_tabix_path,
    header_names = c("chr", "pos", "val"),
    overwrite = TRUE
  ))
  tabix_named_cols <- DBI::dbGetQuery(con, "PRAGMA table_info('tabix_named')")
  expect_equal(tabix_named_cols$name, c("chr", "pos", "val"))

  expect_silent(rduckhts_tabix(
    con,
    "tabix_auto",
    meta_tabix_path,
    auto_detect = TRUE,
    overwrite = TRUE
  ))
  tabix_auto_type <- DBI::dbGetQuery(
    con,
    "SELECT typeof(column1) AS t FROM tabix_auto LIMIT 1"
  )$t[1]
  expect_equal(tabix_auto_type, "BIGINT")

  expect_silent(rduckhts_tabix(
    con,
    "tabix_types",
    meta_tabix_path,
    column_types = c("VARCHAR", "BIGINT", "VARCHAR"),
    overwrite = TRUE
  ))
  tabix_types_type <- DBI::dbGetQuery(
    con,
    "SELECT typeof(column1) AS t FROM tabix_types LIMIT 1"
  )$t[1]
  expect_equal(tabix_types_type, "BIGINT")

  expect_true(DBI::dbExistsTable(con, "vep_variants"))
  if (DBI::dbExistsTable(con, "vep_variants")) {
    expect_silent(DBI::dbGetQuery(
      con,
      "SELECT VEP_Allele FROM vep_variants LIMIT 1"
    ))
  }

  encoded_roundtrip <- DBI::dbGetQuery(
    con,
    paste(
      "SELECT seq_decode_4bit(seq_encode_4bit('ACGTRYSWKMBDHVN')) AS seq,",
      "       seq_encode_4bit('ACGU')::VARCHAR AS u_encode,",
      "       seq_decode_4bit([1::UTINYINT, 0::UTINYINT]) AS eq_decode,",
      "       seq_encode_4bit('ACG!')::VARCHAR AS unknown_encode,",
      "       seq_decode_4bit([1::UTINYINT, 16::UTINYINT]) IS NULL AS invalid_decode"
    )
  )
  expect_equal(encoded_roundtrip$seq[1], "ACGTRYSWKMBDHVN")
  # U (RNA) normalizes to T (code 8)
  expect_equal(encoded_roundtrip$u_encode[1], "[1, 2, 4, 8]")
  # code 0 ('=') is valid nt16
  expect_equal(encoded_roundtrip$eq_decode[1], "A=")
  # unknown chars map to N (15) — permissive like htslib
  expect_equal(encoded_roundtrip$unknown_encode[1], "[1, 2, 4, 15]")
  # out-of-range code (>15) yields NULL
  expect_true(encoded_roundtrip$invalid_decode[1])

  # Test overwrite parameter validation
  if (DBI::dbExistsTable(con, "test_table")) {
    DBI::dbRemoveTable(con, "test_table")
  }

  expect_silent(rduckhts_bcf(con, "test_table", bcf_path, overwrite = TRUE))
  expect_error(rduckhts_bcf(con, "test_table", bcf_path, overwrite = FALSE))

  message("Table creation pattern tests passed!")
}

# Run the test
test_table_creation()

message("Integration tests completed!")