python
seq="MMSRTPSAEPPSDPSSGPPSDPSSGPSSGPPAAQRPLRVLTARHRHEAHRAATPLELFFDLCFVVAIAQAGVQLVHSVAESHAGEGVLDYAMVFFAIWWAWMNFTWFASAYDNDDVPYRLVTLVQIAGVLVLAAGVSRAFEEHDWTAVVIGYVIMRVALSAQWLRAARCAPPPERAMALKYAGGVLVCQVGWLGLLLLPEGSRGWLFLVMVVLELSVPVLAERRHRSSWHPHHIAERYGLLTIIVLGETIAAATVAVKSGIDENDALGELLPIAAGGLLIVFAAWWVYFVVPVHGHLRSSRQAFLWGYGHYLVFASAAAIGAGLEVAVEQAVGEAHISTLAASASVTLPTALFLLTVWALHARHFKVGVAEQLVLPTAALLVVCCTFVGDLAVLAAGLVLAAAVATGVMLTARRTAPGRARTTPAS"
len(seq)


426

python
kd = {'A':1.8,'C':2.5,'D':-3.5,'E':-3.5,'F':2.8,'G':-0.4,'H':-3.2,'I':4.5,'K':-3.9,'L':3.8,'M':1.9,'N':-3.5,'P':-1.6,'Q':-3.5,'R':-4.5,'S':-0.8,'T':-0.7,'V':4.2,'W':-0.9,'Y':-1.3}
import numpy as np
w=19
vals=np.array([kd[a] for a in seq])
window=np.convolve(vals, np.ones(w)/w, mode='same')
threshold=1.6
# identify segments above threshold for >=18
above=window>threshold
segments=[]
i=0
while i<len(seq):
    if above[i]:
        j=i
        while j<len(seq) and above[j]:
            j+=1
        if j-i>=18:
            segments.append((i+1,j))  # 1-based inclusive start, exclusive end
        i=j
    else:
       

([(374, 404)], 1)

python
def predict(th=1.0, minlen=17):
    above=window>th
    seg=[]
    i=0
    while i<len(seq):
        if above[i]:
            j=i
            while j<len(seq) and above[j]:
                j+=1
            if j-i>=minlen:
                seg.append((i+1,j))
            i=j
        else:
            i+=1
    return seg
for th in [0.8,1.0,1.2,1.4]:
    print(th, predict(th))

python
def scan(avg_th=1.6, win=20):
    scores=[]
    for i in range(len(seq)-win+1):
        scores.append((i+1, np.mean([kd[a] for a in seq[i:i+win]])))
    high=[s for s in scores if s[1]>avg_th]
    # merge consecutive starts into segments
    seg=[]
    if not high: return seg
    starts=[h[0] for h in high]
    i=0
    while i<len(starts):
        s=starts[i]
        e=s
        while i+1<len(starts) and starts[i+1]==starts[i]+1:
            i+=1
            e=starts[i]
        seg.append((s, e+win-1))
        i+=1
    return seg
scan(1.4,20), len(scan(1.4,20))


([(51, 80),
  (115, 141),
  (177, 198),
  (182, 201),
  (193, 215),
  (201, 223),
  (237, 261),
  (264, 297),
  (309, 328),
  (312, 332),
  (336, 360),
  (343, 362),
  (365, 413)],
 13)

python
segs=scan(1.4,20)
# merge overlaps
segs_sorted=sorted(segs)
merged=[]
for s,e in segs_sorted:
    if not merged or s>merged[-1][1]+1:
        merged.append([s,e])
    else:
        merged[-1][1]=max(merged[-1][1],e)
merged


[[51, 80],
 [115, 141],
 [177, 223],
 [237, 261],
 [264, 297],
 [309, 332],
 [336, 362],
 [365, 413]]

python
import matplotlib.pyplot as plt, matplotlib.patches as patches
tms=merged
N_inside=True
# build loops between TMs and termini
# orientation alternates: if N inside, first TM goes inside->out, so top is out, bottom is in? In figure, top is out, bottom in.
# We'll map: for each TM, bottom side corresponds to cytosol (in) if starts from N terminus inside and alternates.
sides=[]
current_in=N_inside
for idx,(s,e) in enumerate(tms,1):
    # entering TM from current_in side, exiting opposite
    enter='in' if current_in else 'out'
    exit_='out' if current_in else 'in'
    sides.append((enter, exit_))
    current_in = not current_in
sides

[('in', 'out'),
 ('out', 'in'),
 ('in', 'out'),
 ('out', 'in'),
 ('in', 'out'),
 ('out', 'in'),
 ('in', 'out'),
 ('out', 'in')]

python
import matplotlib.pyplot as plt
from matplotlib.patches import FancyBboxPatch, Circle
import numpy as np, textwrap, math, os

seq="MMSRTPSAEPPSDPSSGPPSDPSSGPSSGPPAAQRPLRVLTARHRHEAHRAATPLELFFDLCFVVAIAQAGVQLVHSVAESHAGEGVLDYAMVFFAIWWAWMNFTWFASAYDNDDVPYRLVTLVQIAGVLVLAAGVSRAFEEHDWTAVVIGYVIMRVALSAQWLRAARCAPPPERAMALKYAGGVLVCQVGWLGLLLLPEGSRGWLFLVMVVLELSVPVLAERRHRSSWHPHHIAERYGLLTIIVLGETIAAATVAVKSGIDENDALGELLPIAAGGLLIVFAAWWVYFVVPVHGHLRSSRQAFLWGYGHYLVFASAAAIGAGLEVAVEQAVGEAHISTLAASASVTLPTALFLLTVWALHARHFKVGVAEQLVLPTAALLVVCCTFVGDLAVLAAGLVLAAAVATGVMLTARRTAPGRARTTPAS"
tms=[[51,80],[115,141],[177,223],[237,261],[264,297],[309,332],[336,362],[365,413]]

# determine inside/outside using positive-inside assumption (N-terminus inside) and alternating topology
N_inside=True
sides=[]
current_in=N_inside
for (s,e) in tms:
    enter='in' if current_in else 'out'
    exit_='out' if current_in else 'in'
    sides.append((enter,exit_))
    current_in = not current_in

# residues to mark (charged)
pos=set("KR")
neg=set("DE")
# Cys positions
cys=[i+1 for i,a in enumerate(seq) if a=="C"]

fig,ax=plt.subplots(figsize=(12,4.8))
ax.set_xlim(0, 100)
ax.set_ylim(0, 40)
ax.axis('off')

# membrane lines
ax.plot([2,98],[30,30],color='black',lw=1)
ax.plot([2,98],[10,10],color='black',lw=1)
ax.text(50,32,"out",ha='center',va='bottom',fontsize=12,fontweight='bold')
ax.text(50,8,"in",ha='center',va='top',fontsize=12,fontweight='bold')

# draw helices
n=len(tms)
x_positions=np.linspace(8,92,n)
helix_w=6
helix_h=18
for i,(x,(s,e)) in enumerate(zip(x_positions,tms),1):
    rect=FancyBboxPatch((x-helix_w/2,11),helix_w,helix_h,
                        boxstyle="round,pad=0.02,rounding_size=1.5",
                        linewidth=1,edgecolor='black',facecolor='white')
    ax.add_patch(rect)
    ax.text(x,20,str(i),ha='center',va='center',fontsize=12,fontweight='bold')
    # residue numbers at top/bottom
    ax.text(x,29.5,str(e),ha='center',va='top',fontsize=9)
    ax.text(x,10.5,str(s),ha='center',va='bottom',fontsize=9)

# loops with beads
def loop_points(x1,x2,y,bulge=4,npts=8):
    xs=np.linspace(x1,x2,npts)
    mid=(x1+x2)/2
    ys=y + bulge*np.sin(np.linspace(0,math.pi,npts))
    return xs,ys

# N-terminus and C-terminus positions
ax.text(4,9,"N",fontsize=12,fontweight='bold')
ax.text(96,9,"C",fontsize=12,fontweight='bold')

# connect N to first helix bottom
ax.plot([4,x_positions[0]-helix_w/2],[9,11],color='black',lw=1)

# loops between helices
for i in range(n-1):
    x1=x_positions[i]+helix_w/2
    x2=x_positions[i+1]-helix_w/2
    # decide loop side based on exit of helix i
    side=sides[i][1]
    y=30 if side=='out' else 10
    bulge=4 if side=='out' else -4
    xs,ys=loop_points(x1,x2,y,bulge=bulge,npts=10)
    ax.plot(xs,ys,color='black',lw=1)
    # beads
    for x,yb in zip(xs[1:-1],ys[1:-1]):
        # map to sequence region between helices
        start=tms[i][1]+1
        end=tms[i+1][0]-1
        if end<start:
            aa=None
        else:
            # pick representative residue proportional
            frac=(x-xs[1])/(xs[-2]-xs[1]) if xs[-2]!=xs[1] else 0
            posi=int(start + frac*(end-start))
            aa=seq[posi-1]
        if aa in pos:
            fc='gray'
        elif aa in neg:
            fc='white'
        else:
            fc='black'
        ax.add_patch(Circle((x,yb),0.5,edgecolor='black',facecolor=fc,lw=0.5))

# connect last helix to C
ax.plot([x_positions[-1]+helix_w/2,96],[11,9],color='black',lw=1)

# highlight cysteines
for c in cys:
    # if in TM label near appropriate helix
    for idx,(s,e) in enumerate(tms,1):
        if s<=c<=e:
            x=x_positions[idx-1]
            ax.text(x+3,20,f"Cys{c}",fontsize=8,rotation=90,va='center',ha='left')
            break

title="Predicted topology model (hydropathy scan + positive-inside rule)\n8 transmembrane helices"
ax.set_title(title,fontsize=12)
outpath="/mnt/data/topology_model.png"
plt.tight_layout()
plt.savefig(outpath,dpi=200)
outpath


'/mnt/data/topology_model.png'