diff --git a/Project.toml b/Project.toml index af02b02..d3ee795 100644 --- a/Project.toml +++ b/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" diff --git a/docs/src/documentation.md b/docs/src/documentation.md index 210dded..05a03bf 100644 --- a/docs/src/documentation.md +++ b/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: diff --git a/docs/src/examples.md b/docs/src/examples.md index 326bdd2..167c644 100644 --- a/docs/src/examples.md +++ b/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)) ``` diff --git a/src/BioStructures.jl b/src/BioStructures.jl index 43329d2..b2976ae 100644 --- a/src/BioStructures.jl +++ b/src/BioStructures.jl @@ -14,6 +14,7 @@ using BioSymbols using CodecZlib using Downloads using Format +import JSON using PrecompileTools using RecipesBase diff --git a/src/download.jl b/src/download.jl index c6418d6..b8557d0 100644 --- a/src/download.jl +++ b/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; ) + downloadalphafold(ids::AbstractArray{<:AbstractString, 1}; ) + downloadalphafold(f::Function, args...; ) + +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; ) + +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 diff --git a/test/runtests.jl b/test/runtests.jl index 254b432..8271e9a 100644 --- a/test/runtests.jl +++ b/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