Add AlphaFold DB download functions

Project.toml

@@ -9,6 +9,7 @@ BioSymbols = "3c28c6f8-a34d-59c4-9654-267d177fcfa9"
 CodecZlib = "944b1d66-785c-5afd-91f1-9de20f533193"
 Downloads = "f43a241f-c20a-4ad4-852c-f6b1247861c6"
 Format = "1fa38f19-a742-5d3f-a2b9-30dd87b9d5f8"
+JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
@@ -44,6 +45,7 @@ DataFrames = "1"
 Downloads = "1"
 Format = "1.3"
 Graphs = "1"
+JSON = "0.21, 1"
 LinearAlgebra = "1.9"
 MMTF = "1"
 MetaGraphs = "0.7, 0.8, 0.9"

docs/src/documentation.md

@@ -683,6 +683,29 @@ end
 Note that some PDB entries, e.g. large viral assemblies, are not available as PDB format files.
 In this case download the mmCIF file instead.
 
+## Downloading AlphaFold models
+
+Predicted structures from the [AlphaFold Protein Structure Database](https://alphafold.ebi.ac.uk) (AlphaFold DB) are identified by UniProt accession.
+To download the latest version of a model:
+
+```julia
+downloadalphafold("P04637", dir="path/to/model/directory")
+```
+
+The file keeps its AlphaFold DB name, e.g. `AF-P04637-F1-model_v6.pdb`, which records the model version.
+[`downloadalphafold`](@ref) accepts the same `format` argument (PDBFormat or MMCIFFormat), lists of IDs and function arguments as [`downloadpdb`](@ref).
+When given a list, IDs with no AlphaFold DB entry are skipped with a warning.
+
+To download and read a model in one step:
+
+```julia-repl
+julia> struc = retrievealphafold("P04637")
+MolecularStructure AF-P04637-F1-model_v6.cif with 1 models, 1 chains (A), 393 residues, 3060 atoms
+```
+
+The per-residue confidence score (pLDDT) is stored as the temperature factor of each atom and can be accessed with [`tempfactor`](@ref).
+[`alphafoldentry`](@ref) returns the AlphaFold DB metadata for an entry, such as the model version, the UniProt sequence and the mean pLDDT.
+
 ## Reading PDB files
 
 To parse an existing PDB file into a Structure-Model-Chain-Residue-Atom framework:

docs/src/examples.md

@@ -136,8 +136,7 @@ CSV.write("1ALW.csv", df) # CSV file writing
 
 ```julia
 using DSSP_jll
-fp = download("https://alphafold.ebi.ac.uk/files/AF-P24941-F1-model_v4.pdb")
-struc = read(fp, PDBFormat; run_dssp=true)
+struc = retrievealphafold("P24941"; run_dssp=true)
 sscodes = sscode.(collectresidues(struc))
 println(join(sscodes))
 ```

src/BioStructures.jl

@@ -14,6 +14,7 @@ using BioSymbols
 using CodecZlib
 using Downloads
 using Format
+import JSON
 using PrecompileTools
 using RecipesBase
 

src/download.jl

@@ -7,9 +7,13 @@ export
     downloadentirepdb,
     updatelocalpdb,
     downloadallobsoletepdb,
-    retrievepdb
+    retrievepdb,
+    alphafoldentry,
+    downloadalphafold,
+    retrievealphafold
 
 const pdb_download_prefix = "https://files.wwpdb.org/pub/pdb"
+const alphafold_api_prefix = "https://alphafold.ebi.ac.uk/api/prediction"
 
 """
     pdbentrylist()
@@ -412,6 +416,157 @@ function retrievepdb(pdbid::AbstractString;
     read(pdbpath, format; structure_name=structure_name, kwargs...)
 end
 
+# The AlphaFold DB file server rejects requests whose User-Agent begins with
+# "curl/", which is the Downloads.jl default
+alphafold_headers() = ["User-Agent" => "BioStructures.jl/$(pkgversion(@__MODULE__))"]
+
+"""
+    alphafoldentry(id::AbstractString)
+
+Obtain the metadata for a predicted structure in the AlphaFold Protein Structure
+Database (AlphaFold DB).
+
+`id` is a UniProt accession, e.g. `"P04637"`, optionally with an isoform
+suffix, e.g. `"P04637-2"`, or an AlphaFold DB entry ID, e.g. `"AF-P04637-F1"`.
+Returns a dictionary with fields that include `"entryId"`, `"latestVersion"`,
+`"uniprotSequence"`, `"globalMetricValue"` (the mean pLDDT), `"pdbUrl"` and
+`"cifUrl"`; see https://alphafold.ebi.ac.uk/api-docs for the full list.
+Throws a `RequestError` if AlphaFold DB has no entry for `id`.
+Requires an internet connection.
+"""
+function alphafoldentry(id::AbstractString)
+    if !occursin(r"^[a-zA-Z0-9-]+$", id)
+        throw(ArgumentError("Not a valid UniProt accession or AlphaFold DB ID: \"$id\""))
+    end
+    io = IOBuffer()
+    Downloads.download("$alphafold_api_prefix/$id", io; headers=alphafold_headers())
+    entries = JSON.parse(String(take!(io)))
+    # A query for a canonical accession also returns entries for its isoforms
+    entryid = uppercase(startswith(uppercase(id), "AF-") ? id : "AF-$id-F1")
+    i = findfirst(entry -> entry["entryId"] == entryid, entries)
+    if i === nothing
+        error("AlphaFold DB returned no entry with ID $entryid; returned IDs: ",
+              join((entry["entryId"] for entry in entries), ", "))
+    end
+    return entries[i]
+end
+
+"""
+    downloadalphafold(id::AbstractString; <keyword arguments>)
+    downloadalphafold(ids::AbstractArray{<:AbstractString, 1}; <keyword arguments>)
+    downloadalphafold(f::Function, args...; <keyword arguments>)
+
+Download the latest version of a predicted structure from the AlphaFold Protein
+Structure Database (AlphaFold DB).
+
+When given an `AbstractString`, a UniProt accession or AlphaFold DB entry ID as
+described for [`alphafoldentry`](@ref), downloads the model file and returns the
+path to the file.
+The file keeps its AlphaFold DB name, e.g. `"AF-P04637-F1-model_v6.pdb"`, which
+records the model version.
+AlphaFold DB serves only the latest version of each model.
+When given an `Array{<:AbstractString, 1}`, downloads the model files in the
+array and returns an array of the paths to the files; IDs with no AlphaFold DB
+entry are skipped with a warning.
+When given a function as the first argument, downloads the file(s) to a
+temporary directory, runs the function with the downloaded filepath(s) as an
+argument, then removes the file(s).
+Requires an internet connection.
+
+# Arguments
+- `dir::AbstractString=pwd()`: the directory to which the file is downloaded;
+    defaults to the current working directory.
+- `format::Type=PDBFormat`: the format of the file; options are PDBFormat and
+    MMCIFFormat.
+- `overwrite::Bool=false`: if set `true`, overwrites the file if it exists in
+    `dir`; by default skips downloading the file if it exists.
+"""
+function downloadalphafold(id::AbstractString;
+                dir::AbstractString=pwd(),
+                format::Type{<:Union{PDBFormat, MMCIFFormat}}=PDBFormat,
+                overwrite::Bool=false)
+    entry = alphafoldentry(id)
+    url = entry[format == PDBFormat ? "pdbUrl" : "cifUrl"]
+    if !isdir(dir)
+        @info "Creating directory: $dir"
+        mkpath(dir)
+    end
+    path = joinpath(dir, split(url, "/")[end])
+    if isfile(path) && !overwrite
+        @info "File exists: $(entry["entryId"])"
+    else
+        @info "Downloading file from AlphaFold DB: $(entry["entryId"])"
+        # Download to a temporary file so that an interrupted download leaves nothing at `path`
+        temppath = tempname(dir)
+        try
+            Downloads.download(url, temppath; headers=alphafold_headers())
+            mv(temppath, path; force=true)
+        finally
+            rm(temppath, force=true)
+        end
+    end
+    return path
+end
+
+function downloadalphafold(ids::AbstractArray{<:AbstractString, 1}; kwargs...)
+    paths = String[]
+    missinglist = String[]
+    for id in ids
+        try
+            push!(paths, downloadalphafold(id; kwargs...))
+        catch err
+            if err isa RequestError && err.response.status == 404
+                push!(missinglist, id)
+            else
+                rethrow()
+            end
+        end
+    end
+    if !isempty(missinglist)
+        @warn "$(length(missinglist)) ID(s) have no AlphaFold DB entry: $missinglist"
+    end
+    return paths
+end
+
+function downloadalphafold(f::Function, args...; kwargs...)
+    mktempdir() do dir
+        f(downloadalphafold(args...; dir, kwargs...))
+    end
+end
+
+"""
+    retrievealphafold(id::AbstractString; <keyword arguments>)
+
+Download and read the latest version of a predicted structure from the AlphaFold
+Protein Structure Database (AlphaFold DB), returning a `MolecularStructure`.
+
+`id` is a UniProt accession or AlphaFold DB entry ID as described for
+[`alphafoldentry`](@ref).
+The per-residue confidence score (pLDDT) is stored as the temperature factor
+of each atom, see [`tempfactor`](@ref).
+Requires an internet connection.
+
+# Arguments
+- `dir::AbstractString=pwd()`: the directory to which the file is downloaded;
+    defaults to the current working directory.
+- `format::Type=MMCIFFormat`: the format of the file; options are PDBFormat and
+    MMCIFFormat.
+- `overwrite::Bool=false`: if set `true`, overwrites the file if it exists in
+    `dir`; by default skips downloading the file if it exists.
+- `structure_name::AbstractString`: the name given to the returned
+    `MolecularStructure`; defaults to the file name, e.g.
+    `"AF-P04637-F1-model_v6.cif"`.
+- other keyword arguments, e.g. `run_dssp`, are passed to `read`.
+"""
+function retrievealphafold(id::AbstractString;
+            dir::AbstractString=pwd(),
+            format::Type{<:Union{PDBFormat, MMCIFFormat}}=MMCIFFormat,
+            overwrite::Bool=false,
+            kwargs...)
+    path = downloadalphafold(id; dir, format, overwrite)
+    read(path, format; kwargs...)
+end
+
 const FMT_SPECIMEN_PATH = joinpath(dirname(dirname(@__FILE__)), "BioFmtSpecimens")
 
 # Get the path to BioFmtSpecimens and download it if required

test/runtests.jl

@@ -114,7 +114,64 @@ end
 
 Aqua.test_all(BioStructures; ambiguities=(recursive=false))
 
-# This is the only test set that requires an internet connection
+# These are the only test sets that require an internet connection
+@testset "AlphaFold DB interface" begin
+    entry = alphafoldentry("Q7TQA6")
+    @test entry["entryId"] == "AF-Q7TQA6-F1"
+    @test entry["uniprotAccession"] == "Q7TQA6"
+    v = entry["latestVersion"]
+    # A query for a canonical accession also returns its isoforms
+    @test alphafoldentry("P04637")["entryId"] == "AF-P04637-F1"
+    @test alphafoldentry("P04637-2")["entryId"] == "AF-P04637-2-F1"
+    @test alphafoldentry("af-q7tqa6-f1")["entryId"] == "AF-Q7TQA6-F1"
+    # Invalid ID format
+    @test_throws ArgumentError alphafoldentry("Q7TQA6/31-308")
+    # Valid accession format but no AlphaFold DB entry
+    @test_throws RequestError alphafoldentry("A0A000Z999")
+    # Invalid file format
+    @test_throws TypeError downloadalphafold("Q7TQA6", dir=temp_dir, format=PDBXMLFormat)
+
+    path = downloadalphafold("Q7TQA6", dir=temp_dir)
+    @test path == joinpath(temp_dir, "AF-Q7TQA6-F1-model_v$v.pdb")
+    @test isfile(path) && filesize(path) > 0
+    path = downloadalphafold("Q7TQA6", dir=temp_dir, format=MMCIFFormat)
+    @test path == joinpath(temp_dir, "AF-Q7TQA6-F1-model_v$v.cif")
+    @test isfile(path) && filesize(path) > 0
+    # An existing file is not downloaded again unless overwrite is set
+    write(path, "placeholder")
+    downloadalphafold("Q7TQA6", dir=temp_dir, format=MMCIFFormat)
+    @test read(path, String) == "placeholder"
+    downloadalphafold("Q7TQA6", dir=temp_dir, format=MMCIFFormat, overwrite=true)
+    @test filesize(path) > length("placeholder")
+
+    # Download multiple files, skipping IDs with no entry
+    paths = @test_logs (:warn, r"A0A000Z999") match_mode=:any downloadalphafold(
+                ["Q7TQA6", "A0A000Z999", "P69905"], dir=temp_dir)
+    @test paths == [joinpath(temp_dir, "AF-Q7TQA6-F1-model_v$v.pdb"),
+                    joinpath(temp_dir, "AF-P69905-F1-model_v$(alphafoldentry("P69905")["latestVersion"]).pdb")]
+    # Errors other than a missing entry are not skipped
+    @test_throws ArgumentError downloadalphafold(["Q7TQA6", "Q7TQA6/31-308"], dir=temp_dir)
+
+    # Test function as first argument to downloadalphafold
+    tmppath = Ref("")
+    @test downloadalphafold("Q7TQA6") do fp
+        tmppath[] = fp
+        s = read(fp, PDBFormat)
+        return countresidues(s, standardselector)
+    end == 331
+    @test !isfile(tmppath[])
+    @test downloadalphafold(fps -> length(fps), ["Q7TQA6", "P69905"]) == 2
+
+    # Test retrievealphafold
+    struc = retrievealphafold("Q7TQA6", dir=temp_dir)
+    @test structurename(struc) == "AF-Q7TQA6-F1-model_v$v.cif"
+    @test countresidues(struc, standardselector) == 331
+    @test all(a -> 0 <= tempfactor(a) <= 100, collectatoms(struc))
+    struc = retrievealphafold("Q7TQA6", dir=temp_dir, format=PDBFormat, structure_name="New name")
+    @test structurename(struc) == "New name"
+    @test countresidues(struc, standardselector) == 331
+end
+
 @testset "PDB interface" begin
     @test length(pdbentrylist()) > 100000