
# Launch PyMOL with the "-R" option (pymol -R), 
# open a Python 3 notebook, then execute 
# remote control of PyMOL (RPC) with XML-RPC:
# https://pymolwiki.org/index.php/Jupyter
# https://pymolwiki.org/index.php/RPC

import xmlrpc.client as xmlrpclib
cmd = xmlrpclib.ServerProxy('http://localhost:9123')

from IPython import display


# PDB: 1kskA01, 1vioA01, 1vs5D02, 2istA01, 2janA03, 2k6pA00, 2xzmW01, 3hp7A01, 3kbgA01, 5mmjd02

cmd.load('/home/eva/pymol/1kskA01.pdb')
cmd.load('/home/eva/pymol/1vioA01.pdb')
cmd.load('/home/eva/pymol/1vs5D02.pdb')
cmd.load('/home/eva/pymol/2istA01.pdb')
cmd.load('/home/eva/pymol/2janA03.pdb')
cmd.load('/home/eva/pymol/2k6pA00.pdb')
cmd.load('/home/eva/pymol/2xzmW01.pdb')
cmd.load('/home/eva/pymol/3hp7A01.pdb')
cmd.load('/home/eva/pymol/3kbgA01.pdb')
cmd.load('/home/eva/pymol/5mmjd02.pdb')


structures = ["1vioA01" ,"1vs5D02", "2istA01", "2janA03", "2k6pA00", "2xzmW01", "3hp7A01", "3kbgA01", "5mmjd02"]
reference = "1kskA01"

for struct in structures:
    cmd.super(struct, reference)
    
cmd.orient(reference)


cmd.bg_color('white')


display.Image('/home/eva/pymol/super_10_1.png')


display.Image('/home/eva/pymol/super_10_grid_1.png')


display.Image('/home/eva/pymol/super_10_2.png')


display.Image('/home/eva/pymol/super_10_grid_2.png')


cmd.color('gray80', '(elem C and model 1kskA01)')
cmd.color('yelloworange', '(elem C and model 1vs5D02)')
cmd.do('''
set cartoon_transparency, 0.6
set ray_trace_mode, 1
''')


label sele and name CA, "%s-%s" % (resn, resi)
set label_position, [0.0, 0.0, 5.0]
set label_size, 10


display.Image('/home/eva/pymol/1kskA01_1vs5D02_overall.png')


display.Image('/home/eva/pymol/1kskA01_beta.png')


display.Image('/home/eva/pymol/1vs5D02_beta.png')


from sys import argv
from collections import defaultdict, Counter
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd


def ss_simplify(code):
    if code in ['E','H']:
        return code
    else:
        return 'C'

def get_ss(dssp_fn, output_fn):
    out = open(output_fn, 'w')
    
    with open(dssp_fn, 'r') as fh:
        for line in fh:
            tk = line.strip().split()
            if tk[-1] == '.' or tk[0] == '#':
                continue
            else:
                resname = line[13]
                if resname in ['a','b']:
                    resname = 'C'
                elif resname in ['X','!']:
                    continue
                else:
                    pass
                resid = int(line[5:10])
                sscode = line[16]
                out.write('%d %s %s\n' % (resid,
                                          resname,
                                          ss_simplify(sscode)))
    
    out.close()
    
    
def heatmap(data, row_labels, col_labels, ax=None,
            cbar_kw={}, cbarlabel="", **kwargs):

    if not ax:
        ax = plt.gca()

    # Plot the heatmap
    im = ax.imshow(data, **kwargs)

    # Create colorbar
    cbar = ax.figure.colorbar(im, ax=ax, **cbar_kw)
    cbar.ax.set_ylabel(cbarlabel, rotation=-90, va="bottom")

    # We want to show all ticks...
    ax.set_xticks(np.arange(data.shape[1]))
    ax.set_yticks(np.arange(data.shape[0]))
    # ... and label them with the respective list entries.
    ax.set_xticklabels(col_labels)
    ax.set_yticklabels(row_labels)

    # Let the horizontal axes labeling appear on top.
    ax.tick_params(top=True, bottom=False,
                   labeltop=True, labelbottom=False)

    # Rotate the tick labels and set their alignment.
    plt.setp(ax.get_xticklabels(), rotation=-30, ha="right",
             rotation_mode="anchor")

    # Turn spines off and create white grid.
    for edge, spine in ax.spines.items():
        spine.set_visible(False)

    ax.set_xticks(np.arange(data.shape[1]+1)-.5, minor=True)
    ax.set_yticks(np.arange(data.shape[0]+1)-.5, minor=True)
    ax.grid(which="minor", color="w", linestyle='-', linewidth=3)
    ax.tick_params(which="minor", bottom=False, left=False)

    return im, cbar


def annotate_heatmap(im, data=None, valfmt="{x:.2f}",
                     textcolors=("black", "white"),
                     threshold=None, **textkw):

    if not isinstance(data, (list, np.ndarray)):
        data = im.get_array()

    # Normalize the threshold to the images color range.
    if threshold is not None:
        threshold = im.norm(threshold)
    else:
        threshold = im.norm(data.max())/2.

    # Set default alignment to center, but allow it to be
    # overwritten by textkw.
    kw = dict(horizontalalignment="center",
              verticalalignment="center")
    kw.update(textkw)

    # Get the formatter in case a string is supplied
    if isinstance(valfmt, str):
        valfmt = matplotlib.ticker.StrMethodFormatter(valfmt)

    # Loop over the data and create a `Text` for each "pixel".
    # Change the text's color depending on the data.
    texts = []
    for i in range(data.shape[0]):
        for j in range(data.shape[1]):
            kw.update(color=textcolors[int(im.norm(data[i, j]) > threshold)])
            text = im.axes.text(j, i, valfmt(data[i, j], None), **kw)
            texts.append(text)

    return texts


path = '/home/eva/pymol/'

dssps = ['1kskA01.dssp',
         '1vioA01.dssp',
         '1vs5D02.dssp',
         '2istA01.dssp',
         '2janA03.dssp',
         '2k6pA00.dssp',
         '2xzmW01.dssp',
         '3hp7A01.dssp',
         '3kbgA01.dssp',
         '5mmjd02.dssp']


for dssp in dssps:
    out = (dssp.replace('dssp', 'txt'))
    get_ss(path + dssp, path + out)


d_ss = defaultdict(list)
d_aa = defaultdict(list)
txts = [dssp.replace('dssp', 'txt') for dssp in dssps]

for txt in txts:
    with open(path + txt) as rf:
        for i in rf:
            d_ss[i.split()[1]].append(i.split()[2])
            d_aa[i.split()[2]].append(i.split()[1])


counted_ss = {i:Counter(j) for i, j in d_ss.items()}
counted_aa = {i:Counter(j) for i, j in d_aa.items()}


def ss_propensity(n_ik, n_i, N_k, N):
    P_ik = (n_ik / n_i) / (N_k / N)
    return P_ik


d = defaultdict(dict)

n_is = {i:sum(j.values()) for i, j in counted_ss.items()}
N_ks = {e:sum(k.values()) for e, k in counted_aa.items()}
N = sum(N_ks.values())

for k, v in counted_ss.items():
    for i, j in v.items():
        n_ik = j
        n_i = n_is.get(k)
        N_k = N_ks.get(i)
        prop = ss_propensity(n_ik, n_i, N_k, N)
        #print('Cклонность аминокислоты {} образовывать {} равна {}'.format(k, i, prop))
        d[k].update({i:prop})


df_prop = pd.DataFrame(d)
df_prop


fig, ax = plt.subplots(figsize=(15, 4), dpi=300)

rows = ['C', 'H', 'E']
columns = ['G', 'S', 'H', 'R', 'L', 'D', 'K', 'F', 'I', 'A', 'Q', 'V', 'E', 'N', 'T', 'P', 'Y', 'W', 'M', 'C']

im, cbar = heatmap(df_prop, rows, columns, ax=ax,
                   cmap="GnBu", cbarlabel="Propensity")

texts = annotate_heatmap(im, valfmt="{x:.3f}")

plt.title('Amino acid secondary structure propensity')
fig.tight_layout()
plt.show()

<__array_function__ internals>:5: UserWarning: Warning: converting a masked element to nan.
/home/eva/miniconda/envs/chemenv/lib/python3.8/site-packages/numpy/core/_asarray.py:85: UserWarning: Warning: converting a masked element to nan.
  return array(a, dtype, copy=False, order=order)

	G	S	H	R	L	D	K	F	I	A	Q	V	E	N	T	P	Y	W	M	C
C	1.458537	1.145823	1.682927	0.961672	0.791966	1.217436	1.219512	0.974326	0.461980	0.774681	1.039455	0.545173	1.112443	1.298258	0.868607	1.402439	0.925610	1.121951	1.070953	0.841463
H	0.401163	0.819396	NaN	1.444186	1.415869	1.024245	0.767442	1.013464	1.699042	1.910299	1.274282	0.745824	0.652739	0.687708	0.621155	0.160465	0.722093	0.534884	0.875264	1.203488
E	0.252439	0.751946	NaN	0.649129	1.187948	0.322263	0.585366	1.062901	1.880918	0.721254	0.593974	2.631055	1.026871	0.432753	1.791503	0.673171	1.514634	1.121951	0.917960	1.262195

Задание 1. Знакомство с укладками¶

Задание 2. Работа с разметкой вторичной структуры в ручном режиме¶

Задание 3. Работа с разметкой вторичной структуры в автоматическом режиме¶