プラスミド工学ワークベンチ構築法

本チュートリアルでは、SpliceCraft の核となるアイデアをインタラクティブなノートブック環境内で再現する、Google Colab ネイティブのプラスミド作業台を作成します。ターミナルベースの TUI に依存するのではなく、Biopython、NumPy、Matplotlib を活用して、プラスミドレコードの読み込み、注釈付きゲノム特徴量の正規化、円形および直線状のプラスミドマップの描画、配列統計量の計算、制限酵素切断部位の解析、バーチャル消化（ディジェスト）シミュレーション、オープンリーディングフレーム（ORF）のスキャン、CDS 特徴量のアミノ酸への変換、プライマー設計、そしてプログラムによる配列編集を行います。

外部依存なしでワークフローを確実に実行できるよう、まずは合成されたオフラインのプラスミドから始めます。ただし、必要に応じて NCBI GenBank からの取得やローカルの GenBank ファイルアップロードにも対応しています。

try:

import Bio

except ImportError:

import subprocess, sys

subprocess.run([sys.executable, "-m", "pip", "install", "-q", "biopython"], check=True)

import math, io, textwrap, random

import numpy as np

import matplotlib

import matplotlib.pyplot as plt

from matplotlib.patches import Polygon

from Bio import Entrez, SeqIO

from Bio.Seq import Seq

from Bio.SeqRecord import SeqRecord

from Bio.SeqFeature import SeqFeature, FeatureLocation

from Bio.Restriction import RestrictionBatch, Analysis, CommOnly

try:

from Bio.SeqUtils import MeltingTemp as _mt

except Exception:

_mt = None

EMAIL = "you@example.com"

ACCESSION = "L09137"

USE_NCBI = False

Entrez.email = EMAIL

TYPE_COLORS = {

"CDS": "#d1495b", "gene": "#c1666b", "rep_origin": "#2e86ab",

"promoter": "#e3a008", "terminator": "#8d99ae", "misc_feature": "#6a994e",

"primer_bind": "#9d4edd", "protein_bind": "#457b9d", "source": "#adb5bd",

}

def _color(ftype): return TYPE_COLORS.get(ftype, "#6a994e")

MCS = "GAATTCGAGCTCGGTACCCGGGGATCCTCTAGAGTCGACCTGCAGGCATGCAAGCTT"

def _synthetic_plasmid():

random.seed(19)

body = "".join(random.choice("ACGT") for _ in range(2686))

seq = body[:400] + MCS + body[400:]

rec = SeqRecord(Seq(seq), id="pDEMO", name="pDEMO",

description="Synthetic offline demo plasmid (SpliceCraft-Colab)")

rec.annotations["topology"] = "circular"

rec.annotations["molecule_type"] = "ds-DNA"

def feat(s, e, strand, ftype, label):

f = SeqFeature(FeatureLocation(s, e, strand=strand), type=ftype)

f.qualifiers["label"] = [label]; return f

rec.features += [

feat(400, 400 + len(MCS), 1, "misc_feature", "MCS"),

feat(700, 1561, 1, "CDS", "AmpR"),

feat(1700, 2260, -1, "rep_origin", "ori"),

feat(2350, 2686, 1, "CDS", "lacZ-alpha"),

]

rec.features[1].qualifiers["product"] = ["beta-lactamase (demo)"]

return rec

def load_record():

if USE_NCBI and "@" in EMAIL and not EMAIL.startswith("you@"):

try:

h = Entrez.efetch(db="nucleotide", id=ACCESSION,

rettype="gbwithparts", retmode="text")

rec = SeqIO.read(h, "genbank"); h.close()

print(f"✓ Fetched {ACCESSION} from NCBI: {rec.description}")

return rec

except Exception as e:

print(f"image NCBI fetch failed ({e}); using synthetic demo plasmid.")

else:

print("image USE_NCBI is off (or EMAIL not set) — using synthetic demo plasmid.")

return _synthetic_plasmid()

def upload_gb():

from google.colab import files

up = files.upload()

name = next(iter(up))

return SeqIO.read(io.StringIO(up[name].decode()), "genbank")

def label_of(f):

for k in ("label", "gene", "product", "note"):

if k in f.qualifiers: return f.qualifiers[k][0]

return f.type

def norm_features(rec, skip_source=True, min_len=0):

L = len(rec.seq); out = []

for f in rec.features:

if skip_source and f.type in ("source",): continue

s = int(f.location.start); e = int(f.location.end)

if (e - s) < min_len: continue

out.append(dict(start=s % L, end=e % L, strand=f.location.strand or 1,

type=f.type, label=label_of(f), color=_color(f.type)))

return out

ノートブック環境を設定し、必要に応じて Biopython をインストールします。さらに、ワークフローに必要な科学計算ライブラリ、可視化ライブラリ、配列解析ライブラリをインポートします。

次に、プラスミドの設定情報や機能ごとの色パレット、多重クローニングサイトの定義、そしてインターネット接続がない場合でもチュートリアルが実行できるよう用意したフォールバック用の合成プラスミドを設定します。また、NCBI や GenBank ファイルからレコードを読み込むためのヘルパー関数を作成し、注釈付き配列特徴を描画可能なメタデータに変換する処理も実装しています。

def _ang(bp, L): return math.pi/2 - 2*math.pi*(bp/L)

def _pt(bp, r, L):

a = _ang(bp, L); return r*math.cos(a), r*math.sin(a)

def _arc(s, e, r, L, n=240):

if e < s: e += L

bps = np.linspace(s, e, max(2, int(n*(e-s)/L)+2))

return ([r*math.cos(_ang(b, L)) for b in bps],

[r*math.sin(_ang(b, L)) for b in bps])

def gc_percent(seq):

s = str(seq).upper(); n = len(s) or 1

return 100.0 * (s.count("G") + s.count("C")) / n

def circular_map(rec, title=None, show_gc=True):

L = len(rec.seq); feats = norm_features(rec)

fig, ax = plt.subplots(figsize=(8, 8)); R = 1.0

if show_gc:

w = max(30, L // 120); step = max(1, w // 2); mean = gc_percent(rec.seq)

s = str(rec.seq).upper()

for i in range(0, L, step):

win = s[i:i+w] or s[i:] + s[:(i+w) % L]

dev = (gc_percent(win) - mean) / 100.0

rr = 0.72 + dev * 0.9

x0, y0 = _pt(i, 0.72, L); x1, y1 = _pt(i, rr, L)

ax.plot([x0, x1], [y0, y1],

color="#2e86ab" if dev >= 0 else "#d1495b", lw=1, alpha=0.5, zorder=1)

xs, ys = _arc(0, L, R, L); ax.plot(xs, ys, color="#333", lw=2.5, zorder=2)

step = max(1, round(L/12/100)*100) or max(1, L//12)

for t in range(0, L, step):

x0, y0 = _pt(t, R*1.015, L); x1, y1 = _pt(t, R*1.045, L)

ax.plot([x0, x1], [y0, y1], color="#999", lw=1, zorder=2)

lx, ly = _pt(t, R*1.10, L)

ax.text(lx, ly, f"{t:,}", ha="center", va="center", fontsize=7, color="#777")

for f in feats:

outer = f["strand"] >= 0

rr = R + 0.09 if outer else R - 0.09

xs, ys = _arc(f["start"], f["end"], rr, L)

ax.plot(xs, ys, color=f["color"], lw=10, solid_capstyle="butt", alpha=0.9, zorder=3)

tip = f["end"] if outer else f["start"]

a = _ang(tip, L); tx, ty = rr*math.cos(a), rr*math.sin(a)

d = -1 if outer else 1

tanx, tany = -math.sin(a)*d, math.cos(a)*d

px, py = math.cos(a), math.sin(a)

ln, wd = 0.055, 0.052

ax.add_patch(Polygon([(tx+tanx*ln, ty+tany*ln),

(tx+px*wd, ty+py*wd),

(tx-px*wd, ty-py*wd)],

color=f["color"], zorder=4))

span = (f["end"] - f["start"]) % L

mid = (f["start"] + span/2) % L

lx, ly = _pt(mid, (rr + 0.16) if outer else (rr - 0.16), L)

ax.text(lx, ly, f["label"], ha="center", va="center",

fontsize=8.5, color=f["color"], weight="bold", zorder=5)

circ = "circular" if rec.annotations.get("topology") == "circular" else "linear"

ax.text(0, 0.05, title or rec.name, ha="center", va="center", fontsize=15, weight="bold")

ax.text(0, -0.07, f"{L:,} bp · {gc_percent(rec.seq):.1f}% GC · {circ}",

ha="center", va="center", fontsize=9.5, color="#555")

ax.set_xlim(-1.45, 1.45); ax.set_ylim(-1.45, 1.45)

ax.set_aspect("equal"); ax.axis("off"); plt.tight_layout(); plt.show()

このコードは、プラスミドの配列データを受け取り、その構造を視覚的に表現する円形マップを描画するための Python スクリプトです。まず、角度や座標を計算する補助関数 _ang と _pt が定義されています。これらは塩基対（bp）の位置から極座標上の点を導き出す役割を果たします。

次に、配列の特定の区間を描画するための _arc 関数が用意されています。これは、始点と終点が円周をまたぐ場合でも正しく処理できるよう設計されており、滑らかな弧を描くために十分な数の点を生成します。さらに、GC 含有率（グアニンとシトシンの割合）を計算する gc_percent 関数も実装されています。

メインとなる circular_map 関数は、これらの機能を活用してプラスミドの全体像を描画します。引数として受け取った配列オブジェクトから長さを取得し、特徴量（遺伝子やプロモーターなど）を正規化して読み込みます。描画領域は 8x8 インチの正方形に設定され、半径は 1.0 とされています。

GC 含有率の可視化オプションが有効な場合、この関数は配列全体をスライディングウィンドウで分析します。各ウィンドウの GC 含有率が平均値からどれだけずれているかを計算し、その偏差に応じて線の色と長さを動的に変更します。GC 含有量が高い領域は青系、低い領域は赤系で表示され、視覚的に違いが一目でわかるようになっています。

次に、プラスミド本体の円周を描画します。太めの線で円を描き、その外側に塩基対の位置を示す目盛りを配置します。各目盛りのラベルには、対応する塩基対の数が表示されます。これらの目盛りは、配列の長さに応じて間隔が調整され、読みやすさを確保しています。

次に、特徴量（遺伝子や機能領域）の描画処理が行われます。各特徴量は、その方向性（ストランド）に応じて円周の内側または外側に配置されます。太い色付きの弧で表現され、矢印によって転写方向が示されます。矢印は、特徴量の始点または終点に配置され、その向きも正確に計算されています。

各特徴名のラベルは、対応する領域の中央付近に表示されます。外側に描画される場合は円の外側、内側のものは内側に配置され、太字で強調表示されます。これにより、複雑なプラスミド構造の中でも重要な情報がすぐに把握できるようになっています。

最後に、図のタイトルとメタデータが表示されます。タイトルが指定されていない場合は、配列の名前が使用されます。また、図の下部には、プラスミドの全長（塩基対数）、GC 含有率、そしてトポロジー（環状または直鎖状）が簡潔に記載されています。描画範囲は適切に調整され、軸は非表示にして視覚的なノイズを排除しています。

このスクリプトを実行することで、研究者や技術者はプラスミドの構造を直感的かつ詳細に確認することが可能になります。特に GC 含有率の分布や遺伝子の配置、方向性などを一目で把握できるため、実験計画や解析において非常に有用なツールとなります。

円形プラスミドの可視化は、塩基対の位置を時計回りのリング上の角度座標にマッピングすることで実現しています。GC 含有量の計算やプラスミドバックボーンの描画、目盛りの追加、特徴アークのレンダリングを行い、さらに鎖の向きが明確に分かるように方向矢印を追加します。また、ラベルと中央の配列メタデータも付与するため、Colab 上でも生物学的な情報が豊富でありながら視覚的に読みやすいマップを作成できます。

Copy CodeCopiedUse a different Browser

def linear_map(rec):

L = len(rec.seq); feats = norm_features(rec)

fig, ax = plt.subplots(figsize=(11, 2.6))

ax.hlines(0, 0, L, color="#333", lw=2)

for f in feats:

y = 0.35 if f["strand"] >= 0 else -0.35

s, e = f["start"], f["end"]

segs = [(s, e)] if e >= s else [(s, L), (0, e)]

for a, b in segs:

ax.annotate("", xy=(b if f["strand"] >= 0 else a, y),

xytext=(a if f["strand"] >= 0 else b, y),

arrowprops=dict(arrowstyle="-|>", color=f["color"], lw=7, alpha=0.85))

ax.text((s + ((e - s) % L)/2) % L, y + (0.28 if y > 0 else -0.28), f["label"],

ha="center", va="center", fontsize=8.5, color=f["color"], weight="bold")

ax.set_xlim(-L*0.02, L*1.02); ax.set_ylim(-1, 1)

ax.set_yticks([]); ax.set_xlabel("position (bp)")

ax.set_title(f"{rec.name} — linear view", fontsize=11)

for sp in ("top", "left", "right"): ax.spines[sp].set_visible(False)

plt.tight_layout(); plt.show()

def gc_skew_plot(rec, window=None):

s = str(rec.seq).upper(); L = len(s)

window = window or max(50, L // 100); skew = []; run = 0

for i in range(0, L, window):

w = s[i:i+window]; g, c = w.count("G"), w.count("C")

run += (g - c) / max(1, g + c); skew.append(run)

fig, ax = plt.subplots(figsize=(11, 2.6))

ax.plot(np.arange(len(skew))*window, skew, color="#2e86ab")

ax.axhline(0, color="#aaa", lw=0.8)

ax.set_title(f"Cumulative GC-skew — {rec.name} (min≈replication origin, max≈terminus)",

fontsize=10)

ax.set_xlabel("position (bp)"); ax.set_ylabel("cumulative skew")

plt.tight_layout(); plt.show()

def stats(rec):

seq = str(rec.seq).upper(); L = len(seq)

print(f"── {rec.name} ─────────────────────────────")

print(f"length : {L:,} bp")

print(f"GC content : {gc_percent(seq):.2f} %")

print(f"A/T/G/C : {seq.count('A')}/{seq.count('T')}/{seq.count('G')}/{seq.count('C')}")

print(f"topology : {rec.annotations.get('topology', 'unknown')}")

print(f"features : {len(norm_features(rec))} annotated")

def _is_circular(rec): return rec.annotations.get("topology") == "circular"

def find_sites(rec, enzymes=None):

batch = RestrictionBatch(enzymes) if enzymes else CommOnly

ana = Analysis(batch, rec.seq, linear=not _is_circular(rec))

return {str(k): v for k, v in ana.full().items()}

def unique_cutters(rec, enzymes=None):

hits = find_sites(rec, enzymes)

return sorted([e for e, pos in hits.items() if len(pos) == 1])

def digest_fragments(rec, enzymes):

hits = find_sites(rec, enzymes)

cuts = sorted({p for pos in hits.values() for p in pos})

L = len(rec.seq)

if not cuts: return [L]

if _is_circular(rec):

frags = [cuts[i+1] - cuts[i] for i in range(len(cuts)-1)]

frags.append(L - cuts[-1] + cuts[0])

else:

frags = [cuts[0]] + [cuts[i+1]-cuts[i] for i in range(len(cuts)-1)] + [L - cuts[-1]]

return sorted((f for f in frags if f > 0), reverse=True)

def restriction_report(rec, enzymes=None):

hits = find_sites(rec, enzymes)

cutting = {e: p for e, p in hits.items() if p}

print(f"── Restriction map of {rec.name} "

f"({'circular' if _is_circular(rec) else 'linear'}) ──")

print(f"unique (single) cutters: {', '.join(unique_cutters(rec, enzymes)) or 'none'}")

print("cutting enzymes (site count):")

for e in sorted(cutting, key=lambda x: len(cutting[x])):

print(f" {e:<10} {len(cutting[e])}x at {cutting[e]}")

より直感的に位置情報を確認できるよう、塩基対軸に沿って注釈された特徴を同じ順序で並べた線形プラスミドマップを作成します。次に、複製に関連する構造を示唆する可能性のある方向性を持つヌクレオチドの偏りを把握するため、累積 GC スキューを計算します。さらに、プラスミドの組成や酵素切断部位、ユニークカッター（単一切断酵素）、予想される消化断片などを要約するための配列統計情報と制限酵素解析ユーティリティも追加しています。

Copy CodeCopiedUse a different Browser

def virtual_gel(named_digests):

lanes = [("Ladder", LADDER)] + list(named_digests.items())

fig, ax = plt.subplots(figsize=(1.4*len(lanes)+1, 6))

ax.set_facecolor("#0b1021")

def y(sz): return math.log10(max(sz, 50))

for i, (name, frags) in enumerate(lanes):

for f in frags:

ax.hlines(y(f), i+0.62, i+1.38, color="#e6f1ff",

lw=2 + min(4, f/1500), alpha=0.9)

if name == "Ladder":

ax.text(i+0.5, y(f), f"{f}", color="#9db4d0", ha="right",

va="center", fontsize=7)

ax.text(i+1, y(max(LADDER))+0.15, name, color="#e6f1ff",

ha="center", fontsize=9, rotation=0)

ax.set_xlim(0.2, len(lanes)+0.5); ax.set_ylim(y(80), y(12000))

ax.invert_yaxis(); ax.axis("off")

ax.set_title("Virtual agarose gel", color="#333", fontsize=11)

plt.tight_layout(); plt.show()

def find_orfs(rec, min_aa=50):

seq = rec.seq; L = len(seq); orfs = []

for strand, s in ((+1, seq), (-1, seq.reverse_complement())):

s = str(s)

for frame in range(3):

i = frame

while i < L - 2:

if s[i:i+3] == "ATG":

j = i

while j < L - 2:

if s[j:j+3] in ("TAA", "TAG", "TGA"):

if (j - i)//3 >= min_aa:

orfs.append((strand, frame, i, j+3, (j-i)//3))

i = j; break

j += 3

i += 3

return sorted(orfs, key=lambda o: -o[4])

def translate_feature(rec, label):

for f in rec.features:

if label.lower() in label_of(f).lower() and f.type in ("CDS", "gene"):

sub = f.location.extract(rec.seq)

prot = sub.translate(table=11, to_stop=True)

return str(prot)

return None

def _tm(primer):

if _mt is not None:

try: return float(_mt.Tm_NN(Seq(primer)))

except Exception: pass

p = primer.upper()

return 2*(p.count("A")+p.count("T")) + 4*(p.count("G")+p.count("C"))

def design_primers(rec, start, end, target_tm=60.0, lo=18, hi=30):

seq = str(rec.seq)

def tune(sub, rev=False):

best = None

for n in range(lo, hi+1):

p = sub[-n:] if rev else sub[:n]

if rev: p = str(Seq(p).reverse_complement())

score = abs(_tm(p) - target_tm)

if best is None or score < best[0]: best = (score, p)

return best[1]

fwd = tune(seq[start:start+hi])

rev = tune(seq[end-hi:end], rev=True)

print(f"── Primers to amplify {start}–{end} ({end-start} bp) ──")

print(f" FWD 5'-{fwd}-3' len {len(fwd)} Tm {_tm(fwd):.1f}°C GC {gc_percent(fwd):.0f}%")

print(f" REV 5'-{rev}-3' len {len(rev)} Tm {_tm(rev):.1f}°C GC {gc_percent(rev):.0f}%")

print(f" amplicon: {end-start} bp")

return fwd, rev

仮想アガロースゲルは、DNA ラダー（マーカー）の隣に切断断片のサイズを対数スケールの移動軸上にプロットすることでシミュレーションします。また、プラスミド配列を全 6 つのリーディングフレームでスキャンし、長いオープンリーディングフレーム（ORF）を特定してアミノ酸長に基づいてランク付けを行います。さらに、注釈付き CDS 機能を翻訳したり、選択した領域周辺にプライマーを設計したりすることも可能です。この際、標的となる融解温度（Tm）に合わせてプライマーの長さを調整します。

def edit_sequence(rec, pos, delete=0, insert=""):
   s = str(rec.seq)
   new = s[:pos] + insert + s[pos+delete:]
   out = SeqRecord(Seq(new), id=rec.id, name=rec.name + "_edit",
                   description=rec.description + f" [edit @{pos}:-{delete}+{len(insert)}"])
   out.annotations = dict(rec.annotations)
   shift = len(insert) - delete
   for f in rec.features:
       s0, e0 = int(f.location.start), int(f.location.end)
       if s0 >= pos: s0 += shift
       if e0 >  pos: e0 += shift
       if 0 <= s0 < e0 <= len(new):
           g = SeqFeature(FeatureLocation(s0, e0, strand=f.location.strand), type=f.type)
           g.qualifiers = dict(f.qualifiers); out.features.append(g)
   return out
LIBRARY = {}
def add_to_library(rec): LIBRARY[rec.name] = rec; print(f"＋ added {rec.name} (now {len(LIBRARY)} in library)")
def show_library():
   for name, r in LIBRARY.items():
       print(f"  {name:<16} {len(r.seq):>7,} bp   {gc_percent(r.seq):.1f}% GC")

挿入、削除、置換に対応した純粋なシーケンス編集機能を実装しました。この機能は新しい編集済みの SeqRecord を返すもので、プラスミドのアノテーションを保持しつつ、下流の機能座標を自動的にシフトさせます。これにより、配列変更後も構築体の内部整合性が保たれます。また、ワークフロー中に元の記録と編集後の記録を保存・一覧表示できる、軽量なノートブック用プラスミドライブラリも用意しました。

Copy CodeCopiedUse a different Browser

if __name__ == "__main__":

rec = load_record()

add_to_library(rec)

print("\n[1] 配列統計情報")

stats(rec)

print("\n[2] 環状マップ（SpliceCraft の特徴的なビュー）")

circular_map(rec)

print("\n[3] 線形マップと GC スキュー")

linear_map(rec)

gc_skew_plot(rec)

print("\n[4] 制限酵素解析")

restriction_report(rec, enzymes=["EcoRI", "BamHI", "HindIII", "PstI",

"SalI", "XbaI", "KpnI", "SacI", "SphI"])

print("\n[5] ゲル上のバーチャル消化")

virtual_gel({

"EcoRI": digest_fragments(rec, ["EcoRI"]),

"EcoRI+HindIII": digest_fragments(rec, ["EcoRI", "HindIII"]),

"BamHI+SalI": digest_fragments(rec, ["BamHI", "SalI"]),

})

print("\n[6] 6 フレーム ORF スキャン（長さ上位 5 つ）")

for strand, frame, s, e, aa in find_orfs(rec, min_aa=40)[:5]:

print(f" strand {strand:+d} frame {frame} {s}-{e} {aa} aa")

print("\n[7] CDS 特徴の翻訳")

prot = translate_feature(rec, "AmpR")

if prot: print(f" AmpR タンパク質 ({len(prot)} aa): {prot[:60]}...")

print("\n[8] 増幅用のプライマー設計")

design_primers(rec, 400, 900)

print("\n[9] 配列編集（6 ベースのタグ挿入）と再マップ")

edited = edit_sequence(rec, pos=456, insert="CATCATCAT")

add_to_library(edited)

print("\n[10] ライブラリ表示"); show_library()

print("\n完了。任意の関数を直接呼び出せます。例：circular_map(edited)。")

まずは、プラスミドレコードを読み込んでライブラリに追加し、配列統計情報を表示するフルガイドデモを実行します。その後、順次、円形マップと線形マップの作成、GC スキュープロット、制限酵素解析レポート、バーチャルゲルの生成、ORF スキャン、CDS 翻訳、そしてプライマー設計出力を処理していきます。

最後に、プラスミドの編集を行い、改変されたコンストラクトを保存。さらに、このノートブックが再利用可能なプラスミド作業台として機能する様子も確認します。

結論として、単一のノートブックをコンパクトかつ機能的なプラスミド工学ワークスペースへと変えることで、本チュートリアルは完結しました。配列の読み込みや特徴量抽出から始まり、グラフィカルなプラスミド可視化、GC スキュー解析、制限酵素マッピング、ゲルシミュレーション、ORF 発見、プライマー設計、そして編集可能な配列操作へとステップを進めました。

各処理を再利用可能な Python 関数として実装したおかげで、このワークフローはモジュール化され、実際の GenBank レコードの検証やプラスミド配列の修正、改変コンストラクトの比較、さらにはクローニングや注釈付けのための追加ユーティリティをノートブックに拡張することも可能になりました。全体として、従来のターミナル型プラスミドツールの機能を、再現性が高く視覚的な、そしてノートブックで扱いやすいバイオインフォマティクスワークフローへと変換する方法を示すことができました。