1
$\begingroup$

Does anyone know of a way to programmatically get panther protein class terms and IDs for a gene list, in a way that the classification structure is preserved (PARENT CLASS + CHILD CLASS)?

i.e. For a given gene, I would want to get the list of classes, in away that I can then use the top classification for downstream analyses.

GENEX, PC00244 zinc finger transcription factor > PC00218 DNA-binding transcription factor > PC00264 gene-specific transcriptional regulator

$\endgroup$

1 Answer 1

1
$\begingroup$

One approach is to get the gene's sequence (DNA/AA) and use the interproscan API (e.g., https://www.ebi.ac.uk/Tools/common/tools/help/index.html?tool=iprscan5) with the httr2 library in R to POST gene sequence(s) and GET the resulting Panther annotation. Make sure to include Panther in the appl parameter for the post request.

Here's a function that does this for a single sequence.

submit_fasta <- function(
    fasta,
    path_outfolder,
    email,
    stype = "n",
    appl = c("Panther"),
    title = 'stackoverflow'
    ) {
  stopifnot("stype must be 'n' for nucleotide seqs or 'p' for protein seqs" = stype %in% c('n', 'p'))
  url_post <- "https://www.ebi.ac.uk/Tools/services/rest/iprscan5/run"
  chr_seq <- as.character(fasta)
  req <- httr2::request(url_post)
  req <- req |> httr2::req_body_form(
    "email" = email,
    "title" = title,
    "goterms" = "false",
    "pathways" = "false",
    "stype" = stype,
    "appl" = paste0(appl, collapse = ","),
    "sequence" = chr_seq
  )
  req <- req |> httr2::req_headers("Accept" = "text/plain", "Content-type" = "application/x-www-form-urlencoded")
  resp_post <- req |> httr2::req_perform()
  chr_retrieval_code <- resp_post$body |> rawToChar()
  url_get <- "https://www.ebi.ac.uk/Tools/services/rest/iprscan5/result"
  req <- httr2::request(url_get)
  req <- req |> httr2::req_url_path_append(chr_retrieval_code)
  req <- req |> httr2::req_url_path_append("tsv")
  req <- req |> httr2::req_retry(
    max_tries = 15L, backoff = ~ 30L,
    is_transient = function(resp) {httr2::resp_status(resp) == 400L}
  )
  req <- req |> httr2::req_headers(Accept = NULL)
  req |> httr2::req_dry_run(quiet = FALSE)
  resp_get <- req |> httr2::req_perform()
  cat('status code:', resp_get$status_code, '\n')
  result <- resp_get$body |> rawToChar()
  result <- read.table(text = result, sep = '\t')
  return(result)
}

Importantly, the fasta input parameter is expected to be XString object from the Biostrings package. However, the code could easily be adjusted. The sequence eventually needs to be coerced to a string for the POST request.

#fasta <- Biostrings::readDNAStringSet('seq.fa')
#resp <- submit_fasta(fasta = fasta_single, path_outfolder = 'x', email = '[email protected]', title = 'stackoverflow')

edit: set the is_transient param to handle the actual transient error case where the data is not yet available (response code = 400). Also, it's important to note there's a "status" endpoint https://www.ebi.ac.uk/Tools/common/tools/help/index.html?tool=iprscan5#/Status) which is probably much more appropriate than pinging the data endpoint until either success or reaching max tries. However, I hope this answer suffices!

$\endgroup$
2
  • $\begingroup$ Wow, thank you for your comprehensive answer and welcome to the site. $\endgroup$
    – M__
    Nov 21, 2023 at 13:29
  • 1
    $\begingroup$ @M__ Thanks for the welcome reply! $\endgroup$
    – Jacob Krol
    Nov 23, 2023 at 0:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.