Plannotate Plasmid Annotation

# Install via pip (requires BLAST+ on PATH)
pip install plannotate

# Install via conda (recommended — handles BLAST+ automatically)
conda install -c conda-forge -c bioconda plannotate

# Verify installation
plannotate --help
python -c "import plannotate; print('plannotate OK')"

from plannotate import annotate, write_genbank, create_bokeh_chart
from Bio import SeqIO

# Load plasmid from FASTA
record = next(SeqIO.parse("plasmid.fasta", "fasta"))
sequence = str(record.seq)

# Annotate (circular, against Addgene database)
results = annotate(sequence, linear=False, db="addgene")
print(f"Found {len(results)} features")
print(results[["Feature", "Feature_type", "pct_identity", "pct_query_cov"]].to_string())

# Export GenBank file
write_genbank(sequence, results, output_file="plasmid_annotated.gb")

# Generate interactive HTML map
create_bokeh_chart(sequence, results, output_file="plasmid_map.html")
print("Outputs: plasmid_annotated.gb, plasmid_map.html")

from Bio import SeqIO
import os

# Option A: Load from FASTA
def load_fasta(path):
    record = next(SeqIO.parse(path, "fasta"))
    seq = str(record.seq).upper()
    return seq, record.id

# Option B: Load from GenBank (strip annotations, keep sequence)
def load_genbank(path):
    record = next(SeqIO.parse(path, "genbank"))
    seq = str(record.seq).upper()
    return seq, record.id

# Option C: Raw sequence string
raw_seq = "ATGCGTAAAGGAGAAGAACTTTTCACTGGAGTTGTCCCAATTCTTGTTGAATTAGATGGTGATGTT"

# Validate sequence
def validate_plasmid(seq, name="plasmid"):
    valid_bases = set("ATGCNRYSWKMBDHV")
    invalid = set(seq.upper()) - valid_bases
    if invalid:
        raise ValueError(f"Invalid bases in {name}: {invalid}")
    if len(seq) < 100:
        raise ValueError(f"Sequence too short ({len(seq)} bp); minimum 100 bp")
    gc = (seq.count("G") + seq.count("C")) / len(seq) * 100
    print(f"{name}: {len(seq):,} bp, GC={gc:.1f}%")
    return seq

seq, plasmid_id = load_fasta("plasmid.fasta")
validate_plasmid(seq, plasmid_id)

from plannotate import annotate

# Annotate circular plasmid against the Addgene database (most comprehensive for common vectors)
results = annotate(
    seq,
    linear=False,       # False = circular plasmid (default)
    db="addgene",       # Database: "addgene", "fpbase", or "snapgene"
)

print(f"Total features detected: {len(results)}")
print(f"\nColumns: {list(results.columns)}")

# Preview feature table
cols = ["Feature", "Feature_type", "start", "end", "strand", "pct_identity", "pct_query_cov"]
print(results[cols].sort_values("start").to_string(index=False))

import pandas as pd

# Inspect hit quality distribution
print("Identity percentile summary:")
print(results["pct_identity"].describe().round(1))
print("\nCoverage percentile summary:")
print(results["pct_query_cov"].describe().round(1))

# Separate high- and low-confidence hits
high_conf = results[
    (results["pct_identity"] >= 95) &
    (results["pct_query_cov"] >= 90)
].copy()

low_conf = results[
    (results["pct_identity"] < 95) |
    (results["pct_query_cov"] < 90)
].copy()

print(f"\nHigh-confidence features (identity>=95%, coverage>=90%): {len(high_conf)}")
print(f"Low-confidence / partial features:                       {len(low_conf)}")

if not low_conf.empty:
    print("\nLow-confidence features (review manually):")
    print(low_conf[["Feature", "Feature_type", "pct_identity", "pct_query_cov"]].to_string(index=False))

# Save filtered table
results.to_csv("all_features.csv", index=False)
high_conf.to_csv("high_confidence_features.csv", index=False)
print("\nSaved: all_features.csv, high_confidence_features.csv")

from plannotate import write_genbank

# Write full annotation (all features)
write_genbank(seq, results, output_file="plasmid_annotated.gb")
print("Written: plasmid_annotated.gb")

# Write with high-confidence features only
write_genbank(seq, high_conf, output_file="plasmid_highconf.gb")
print("Written: plasmid_highconf.gb")

# Verify using BioPython
from Bio import SeqIO
record = next(SeqIO.parse("plasmid_annotated.gb", "genbank"))
print(f"\nGenBank verification:")
print(f"  Sequence length: {len(record.seq):,} bp")
print(f"  Features: {len(record.features)}")
for feat in record.features:
    label = feat.qualifiers.get("label", ["(unlabeled)"])[0]
    print(f"  [{feat.type:20s}] {label} @ {feat.location}")

from plannotate import create_bokeh_chart

# Generate interactive circular plasmid map
create_bokeh_chart(
    seq,
    results,
    output_file="plasmid_map.html",
)
print("Interactive map saved: plasmid_map.html")
print("Open in any browser — no server required")

# Tip: open automatically in the default browser
import webbrowser, os
webbrowser.open(f"file://{os.path.abspath('plasmid_map.html')}")

from Bio import SeqIO
from Bio.SeqFeature import FeatureLocation
import pandas as pd

record = next(SeqIO.parse("plasmid_annotated.gb", "genbank"))

# Build a feature DataFrame from the GenBank record
rows = []
for feat in record.features:
    label = feat.qualifiers.get("label", [""])[0]
    note  = feat.qualifiers.get("note",  [""])[0]
    rows.append({
        "type":   feat.type,
        "label":  label,
        "note":   note,
        "start":  int(feat.location.start),
        "end":    int(feat.location.end),
        "strand": feat.location.strand,
        "length": len(feat.location),
    })

feat_df = pd.DataFrame(rows)
print(feat_df.to_string(index=False))

# Example: find antibiotic resistance genes
resistance = feat_df[feat_df["label"].str.contains(
    r"AmpR|KanR|CmR|TetR|SpecR|HygR|ZeoR|BlastR|GentR",
    case=False, na=False, regex=True
)]
print(f"\nAntibiotic resistance markers found: {len(resistance)}")
print(resistance[["label", "start", "end", "length"]].to_string(index=False))

from plannotate import annotate, write_genbank, create_bokeh_chart
from Bio import SeqIO
import pandas as pd
import os

input_dir  = "plasmids/"          # directory of *.fasta files
output_dir = "annotated_results/"
os.makedirs(output_dir, exist_ok=True)

summary_rows = []

for fasta_file in sorted(f for f in os.listdir(input_dir) if f.endswith(".fasta")):
    plasmid_name = fasta_file.replace(".fasta", "")
    fasta_path   = os.path.join(input_dir, fasta_file)

    record = next(SeqIO.parse(fasta_path, "fasta"))
    seq    = str(record.seq).upper()

    print(f"Annotating {plasmid_name} ({len(seq):,} bp)...", end=" ")
    results = annotate(seq, linear=False, db="addgene")
    print(f"{len(results)} features")

    # Save per-plasmid outputs
    write_genbank(
        seq, results,
        output_file=os.path.join(output_dir, f"{plasmid_name}.gb")
    )
    create_bokeh_chart(
        seq, results,
        output_file=os.path.join(output_dir, f"{plasmid_name}.html")
    )
    results.to_csv(os.path.join(output_dir, f"{plasmid_name}_features.csv"), index=False)

    # Accumulate for summary
    results["plasmid"] = plasmid_name
    summary_rows.append(results)

# Consolidated summary table
summary = pd.concat(summary_rows, ignore_index=True)
summary.to_csv(os.path.join(output_dir, "all_plasmids_features.csv"), index=False)
print(f"\nBatch complete. {summary['plasmid'].nunique()} plasmids annotated.")
print(f"Summary table: {output_dir}all_plasmids_features.csv")