Skip to content

Instantly share code, notes, and snippets.

@jbarnoud
Forked from kain88-de/release-0.16.0.ipynb
Last active April 8, 2017 08:18
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save jbarnoud/a87b51205dfa89e9a13547bc3aa3c5d7 to your computer and use it in GitHub Desktop.
Save jbarnoud/a87b51205dfa89e9a13547bc3aa3c5d7 to your computer and use it in GitHub Desktop.
New MDAnalysis features in release 0.16.0
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"run_control": {
"frozen": false,
"read_only": false
}
},
"source": [
"# Attach arbitrary time series to your trajectories"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Attach arbitrary time series to your trajectories\n",
"\n",
"Our GSoC student [@fiona-naughton](https://github.com/fiona-naughton) has implemented an auxillary reader to add arbitrary time series to a universe. The time series are kept in sync with the trajectory so it is possible to iterate through the trajectory and access the auxiliary data corresponding to the current time step."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[ 0. 200.71288 -1552.2849 ..., 128.4072 1386.0378 -2699.3118 ]\n",
"[ 50. -1082.6454 -658.32166 ..., -493.02954 589.8844 -739.2124 ]\n",
"[ 100. -246.27269 146.52911 ..., 484.32501 2332.3767 -1801.6234 ]\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"/home/jon/dev/mdanalysis/package/MDAnalysis/topology/guessers.py:56: UserWarning: Failed to guess the mass for the following atom types: M\n",
" \"\".format(', '.join(misses)))\n"
]
}
],
"source": [
"import MDAnalysis as mda\n",
"from MDAnalysisTests.datafiles import PDB_sub_sol, XTC_sub_sol, XVG_BZ2\n",
"\n",
"# Create your universe as usual\n",
"universe = mda.Universe(PDB_sub_sol, XTC_sub_sol)\n",
"# Attach an auxiliary time serie with the name `forces`\n",
"# In this example, the XVG file contains the force that applies to each atom.\n",
"universe.trajectory.add_auxiliary('forces', XVG_BZ2)\n",
"# Itarete through your trajectory, the time serie is kept in sync\n",
"for time_step in universe.trajectory:\n",
" print(time_step.aux.forces)\n",
"# The first element of each array is the time in picoseconds.\n",
"# The next elements are the other columns of the XVG file."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"[@fiona-naughton](https://github.com/fiona-naughton) worked at offering several convenient way to iterate through your data. Read the documentation or [Fiona’s blog](https://fiona-naughton.github.io/blog/) posts to learn more about the feature.\n",
"\n",
"This feature is still in its beginning and will be expanded in future releases. You can follow the conversation on the [initial issue](https://github.com/MDAnalysis/mdanalysis/issues/785) or on the [pull request](https://github.com/MDAnalysis/mdanalysis/pull/868). So far, only the XVG format used by [gromacs](http://www.gromacs.org/) and [grace](http://plasma-gate.weizmann.ac.il/Grace/) are supported. Open an issue if you need support for other time series formats."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Do a dimension reducion with PCA and Diffusion Maps\n",
"\n",
"[@jdetle](https://github.com/jdetle) has implemented two new dimension reduction algorithms, [Principal Component Analysis](pca) and [Diffusion Maps](dmaps-paper). Both can be found in the analysis submodule. As an example lets look at the first two PCA dimensions of ADK from our test files.\n"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import MDAnalysis as mda\n",
"from MDAnalysis.analysis.pca import PCA\n",
"from MDAnalysisTests.datafiles import PSF, DCD"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Mean Calculation Step 98/98 [100.0%]\n",
"Step 98/98 [100.0%]\n"
]
}
],
"source": [
"u = mda.Universe(PSF, DCD)\n",
"ca = u.select_atoms('protein and name CA')\n",
"\n",
"pca = PCA(u, select='protein and name CA', verbose=True).run()\n",
"reduced_data = pca.transform(ca, n_components=2)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X2cnGV97/HPN5sFNqAuaKSwEBMRQRBNZMsBY3sEraio\nRLSCtUrVI8dz1FOppS6ntj5SttIjtdXTivj4eqGAICsaNYrxqVSsiRsMUXPKMwyPSlaBLLBJfueP\nuSfOTu57ZnZ3Hu575vt+vfaVnZl7Z66d3Du/+7qu3/W7FBGYmZnlzaJuN8DMzCyNA5SZmeWSA5SZ\nmeWSA5SZmeWSA5SZmeWSA5SZmeWSA5RZgUhaLek/JT0kaU2322PWTg5QZg1IulXSdBIU7pX0WUn7\nVT1+sqQfSHpQ0v2Svi/pFTXP8XxJIendC2zOB4CPRcR+ETFRp83fk7RN0t41939W0mNJWx+UdIOk\n8yU9oeqYP5P0b1W3Hy/pWklXStprge03a5oDlFlzXh4R+wHPAUaB9wBIejXwJeDzwCHAgcDfAi+v\n+fkzgQeANyywHU8BttQ7QNJy4A+AAF6RcsiHI+JxwFLgjcDxwLWS9k15rv2B7wC3AadHxGMLabzZ\nXDhAmc1BRJSAbwDPlCTgI8AHI+LiiPhNROyKiO9HxFsqP5N88L8aeBtwuKTReq8h6S2SbpT0gKSr\nJR2c3H8T8FTgq0lvbu+Mp3gDcB3wWcqBMet3eSQifkI5iD2RcrCqbsdS4LvADcCfRsSOeu02azUH\nKLM5kHQo8FJgEjgCOBS4osGPnQY8RLmntY46QUPSScD5wGuAgyj3XC4FiIjDgNtJenMR8WjG07wB\nuCT5OlnSgfUaFxEPAt+m3OuqOAD4HvAj4E0Rsav+r2jWeg5QZs2ZkDQF/BvwfeDvKPc6AO5u8LNn\nApdFxE7gC8AZkgYzjn0d8OmI+GkSgM4FTkiG7RqS9DzKw4CXR8RG4CbgT5r40bsoB6WKQ4GnA58N\nF+y0LnGAMmvOmogYjoinRMT/jIhp4NfJYwdl/VDS4zqRcm8G4CvAPsApGT9yMOVeEwAR8VDyOiNN\ntvNM4FsR8avk9heo02OrMkJ5jqzieuAvgW9IWtXka5u11OJuN8CswLYCdwCvAv4h45jXU74Q/Gp5\nygooB6gzgbQsvLso94CA3fNXTwRKjRojaYjy0OCApHuSu/cGhiU9OyKuz/i5/YAXAudV3x8RH03m\nub4t6fkRcUOjNpi1kntQZvOUDH39BfA3kt6YpGMvkvQ8SRclh50JvB9YWfX1KuClkp6Y8rRfBN4o\naWUSHP4O+HFE3NpEk9YAO4Gjql7rGcAPSckelLS3pGMpB8ptwGdSfscPAx8FrpF0RBNtMGsZByiz\nBYiIK4DTgTdR7v3cC3wI+Iqk4yn3hj4eEfdUfV0N3Ai8NuX5rgH+BriS8tzWYcAZTTbnTOAzEXF7\n9esBHwNeJ6kyYvJXkh6kPHT4eWAj8NyIeDjjd/wgcDHwHUmHNdkWswWT5z/NzCyP3IMyM7NccoAy\nM7NccoAyM7NccoAyM7Nc6ul1UE960pNi+fLl3W6GmZlV2bhx468iYmmj43o6QC1fvpwNGzZ0uxlm\nZlZF0m2Nj/IQn5mZ5ZQDlJmZ5ZIDlJmZ5ZIDlJmZ5ZIDlJmZ5VJPZ/GZ1ZqYLHHBuq3cNTXNwcND\nnHPyEaxZ1exWS2bWSQ5Q1jcmJkuc++XNTM/sBKA0Nc25X94M4CBllkMOUNY3Lli3dXdwqpie2ckF\n67buftw9K7P8cICyvnHX1HTq/ZWelHtWZvniAGV94+DhIUopQWpASu1Zvf+rW9yrMusiZ/FZ3zjn\n5CMYGhyYdd/Q4AA7Mzbt3LZ9htLUNMHvelUTk6WWtGVissTq8fWsGFvL6vH1LXtes17S0zvqjo6O\nhmvx9Za0LDxofv4o7ecvWLc1tWeVZkBiV0Td12mUKVibrAHlQHn+ace4h2Z9QdLGiBhteJwDlBVF\n2gf74IAgYGZXzLpv370W85vpmaaG5tKetxlpQaWZ4LN6fH1qQNx/ySBL9lq8oCFFp9FbEThA4QDV\na7I+2BsREMDIHHo9Dz+6g6npmYbPPTI8xLVjJzVsY3XwafYvrhWB1j0zyyMHKByges2KsbVNf7hn\nafYDey69qpHhod2BbT4BtFlZba8E16zXrg2iZt1W6AAl6dPAy4D7IuKZyX0HAJcBy4FbgddExLZ6\nz+MAVXzVPZtFUmZCw1w0+4HdzGtXemdZt1utdg4MaBhIBdwyfkobW2U2N80GqLymmX8W+Bjw+ar7\nxoDvRMS4pLHk9ru70DbrkNpeTFqASJuDaiRrPVStNatGdvdW0npUacEoMu7PIpjTkGLlPahkFe69\neFHDXt7Bw0NNtsYsX3IZoCLiB5KW19x9KvD85PvPAd/DAaqnpM0DpX34pvUiKj83vGSQhx7ZUTdg\nzecDuxKoqtuXNaRWme9qNJ9V3ZObT6LG9MzOhscPDQ5w4pFLWT2+3okTVji5DFAZDoyIu5Pv7wEO\nTDtI0lnAWQDLli3rUNNsodLq5GXZFbHHkFVtJl1lTqa2NzM0OLA7qM1VdY8KshMiaocQs5IXqttR\nGwCbCbSNjAwPceKRS7lyY8lVMqyQcjkHBZD0oL5WNQc1FRHDVY9vi4j96z2H56Dyr9EEf5q5TPq3\nM+16Lllz82lHM3Ng+y8Z5JGZXZltaDaImnVS0eeg0twr6aCIuFvSQcB93W6QLcx8hrXm2gOq7fW0\nUtqwX1bgmU87Gs2BDQ0O8N6XH123DVnzbc3Ow5l1U5EC1NXAmcB48u9XutscW6i06uK1WrF4tZ3a\nGQBrXweyA1FWG7LmyoaXDHpeynIvlwFK0hcpJ0Q8SdKdwHspB6bLJb0ZuA14TfdaaK3Q6Cq+0kPw\nB2fZfILhOScfkVp946FHdrBtezlxw/NSlle5DFAR8dqMh17Q0YZYW9XLhKtX9cGal9bzSssqrOyL\n5ffb8iSXAcr6Q9rVvUvztF5tz2vF2NrU4zwvZXnjAGUdV52dNrxkkL0XL2q63pwtXFbP1Qt6LW+8\nH5R1VCUbrbLP0rbtMzy6YxcXnr6Sa8dOcnDqgLR9sQYHxMOP7vD+VJYrDlDWUWmZe5X5D+uMNatG\nOP+0YxgZHkKUMyUJmJqeacvmjGbz5QBlHeV1OfmwZtUI146dxC3jp7Bkr8V7VKzwRYPlgeegrK1q\nKyg8YWgwtS6d5z+6xxcNllfuQVnb1M43laamefixHQwu0qzjFlIfzxYu6+LAFw3WbQ5Q1jZp800z\nO4P99lm8e/5jZHjIaeVdlpY04YsGywMP8VlLpBVDzRoimto+w+TfvqjDLbQsaYt5TzxyKRes28rZ\nl21y+r91jQOULVjaVhnnfnmz55sKpF5hWpdCsm7xEJ8tWFbquISHjgrISwEsLxygbMHqDeVVr7fx\nfFMxOKvP8sJDfLZg9UrndGo7CmudrP/PRRIrxtZ6Tso6xj0oWzBngfWWtP9PgJ0RrjRhHeUAZQtW\nWzrHQ3nFVvv/OSDtcYznpKwTFBGNjyqo0dHR2LBhQ7ebYVZoK8bWkvYpIeCW8VM63RzrAZI2RsRo\no+M8B2Xzlrb2yb2m3uPtOaxbPMRn85JWxsjzEr3J23NYtzhA2bx4rUz/8PYc1i0OUDYvXivTX7w9\nh3WDA5TNiytg96+si5DS1LSH/KylHKBsXrz2qX/VuwjxkJ+1kgOUNTQxWWL1+PpZV8de+9S/shby\nVvOQn7WC08ytrrTK1udccT3vu3oLv5me4eDhIS48faUDUx+p3Z4jayWl5yNtoRygrK6sTQcr22h4\nK4b+VF1jcfX4eq+TsrbwEJ/V1cxVsIdz+pvnI61d3IPqc42qQWRVEajl4Zz+5R15rV0coPpY1s6p\nG257gO/+8n7umppmeMkgg4u0x7qXWh7O6W/ekdfawUN8fSyrGsQl192+u4TRtu0zIBgeGtxdRWBw\n0ezq1h7OsWquMmKtUrgelKRbgQeBncCOZiriWrqsYbnavtLMzmDfvRez6b0vAlwk1upzlRFrlcIF\nqMSJEfGrbjei6JqdX4LZHy7eJdfqcfVzaxUP8fWxtOyrPbemK/OHizXLWX3WKkXsQQXwLUkBfCIi\nLqp+UNJZwFkAy5Yt60LziiMr++rKjaVZcwj+cLG5SDuvKufP6vH1Hhq2phVuR11JIxFRkvRk4NvA\nOyLiB2nHekfdPTUzf+Q5Jmu12sw+KF/4uDxWf+rZHXUjopT8e5+kq4DjgNQAZbM1m/7rOSZrtXqZ\nfT7XLEuh5qAk7SvpcZXvgRcBN3S3VcXh9F/rFmf22XwUrQd1IHCVJCi3/QsR8c3uNqk4/CFh3eLM\nPpuPQvWgIuLmiHh28nV0RJzX7TYViTcZtG5xZp/NR9F6UDZH1QkPaWWL/CFhnZCV2bdm1YiTciyT\nA1QPq02K2LZ9hsEBMTw0uHsvJ38YWKekJd+4bp/V4wDVw7L2cqouW2TWTc7us3ocoHpI7VBJVhkj\nJ0VYXjhxx+opVJKEZasMlVSqkJempl22yHLPiTtWjwNUj0gbKgn2rK3npAjLE2f3WT0e4usR9bbO\nGBkecoaU5VK97D4zB6iCq8w7ZVVUHBke4tqxkzraJrO5qM3um5gsuaisAQ5QhZZWgLOah0qsaJx2\nbtU8B1VgafNOFSPDQ64UbYXjepFWzT2oAsuadxJ4WM8KKeucLk1Ns2JsrYf8+ox7UAXmFF3rNfXO\n3cryiXO/vJmJyVLnGmVd4wBVEJWJ4xVja1k9vp6JyZJTdK3npJ3TtTzk1z8coAogbRFuZeL4/NOO\nYWR4COF5Jyu+NatGZp3TWVxpoj94DqoA6k0cXzt2kgOS9ZTqtPPV4+u9j1Qfcw+qAFyvzPqVh7H7\nmwNUATgZwvpV7ZCfh7H7i4f4cqi2KvmJRy7lyo2lWcN8voq0fuFKE/3LPaicSUuIuHJjiVcdO+Kr\nSOt7WQlDTjvvTe5B5UxWQsR3f3m/F99a3/MGh/3FAarLvMmgWfOcMNRfPMTXRd5k0GxunDDUXxyg\nusibDJrNTVra+eCAePjRHbOqrFhvcIDqokabDDohwmy22rTz/ZcMQsDU9IyTJnqQ56C6ZGKyxCKJ\nnbHnVoPeZNAsW22liW3bZ2Y97qSJ3uEeVBdU5p7SgpOH88ya56SJ3uYA1QVZGw0OSB7OM5sDJ030\ntqYDlKQDmvgabmdje0XW1d2uCAcnszlwrb7eNpc5qLuSr3pV8AeAZQtqUQOSXgx8NHmtiyNivJ2v\n1w5Z65181Wc2N5ULuspawuElg0TA2Zdt4oJ1W10GqeDmMsT3i4h4akSsyPoCft2uhgJIGgA+DrwE\nOAp4raSj2vma7eCrPrPWWbNqhGvHTuLC01fyyMwuZ/T1kLkEqBNadMxCHAfcGBE3R8RjwKXAqW1+\nzZZzhWaz1ssqg/Suy6/3GqmCajpARcQjWY9V5p7qHdMiI8AdVbfvTO6rbstZkjZI2nD//fe3uTnz\nU1veyMMQZguXNbe7M8I9qoJqGKAkHSvpvZL2l/Q4ScdLerOkj0haJ6kE3Nr+pjYnIi6KiNGIGF26\ndGm3m7MHV2M2a49m5nAra6SsGJrpQX0C+BpwO7AV+CCwErgROAZYFRGdyt4rAYdW3T4kua8w6lVj\nNrP5S5vbTeM1UsXRTBbfvwPnAD8FlgCfjIjLASSdExH3tbF9tX4CHC5pBeXAdAbwJx18/QXzwkKz\n9qjN6Muq1OJs2eJoGKAi4n9JWhIR2yUdALxH0tnAByiXjeuYiNgh6e3AOspp5p+OiC2dbMNCOcXc\nrH2qyyBVhtO9E3VxNZUkERHbk38fiIi/4Hc9lwMlndjG9qW15esR8fSIOCwizuvka7eCU8zNOsPZ\nssWnSOkCN/3D0kpgHBiKiP/asla1yOjoaGzYsKHbzdiDs/jMrJ9J2hgRo42Oa7qShKSfRsRzqu+L\niE3AiyU9P+sY+53awHTh6SsdmMw6yBeHxTKXUkfPkPSzrAclCXjCwpvUWyp/EJXdciv91Up6OeA/\nELMOqJ2T8t9g/s0lQB3ZxDF7lujuY7V/ELWDqd63xqxz6i3x8N9gPjUdoCLitnY2pBdlbatRzenl\nZp3hJR7F4/2g2qiZE9/p5Wad4b2jiscBqo0anfhOLzfrnLQlHoMD4uFHd7iYbE7NZQ7K5mBissT2\nx3bscX8lUWLEGURmHZW2d9RDj+xganoGcNJEHs07QEn6I+A1wMcjYpOksyLiotY1rbjSVrADDA8N\n8r5XHO2T36xLqitNrB5fz7btM7Med9JEviykB/Um4H9QLn10AOUCskZ2csS+ey/2iW+WE06ayL+F\nzEE9GBFTEfGXwIuA329RmwrPJ75Z/jlpIv8WEqDWVr6JiDHg8wtvTm/wiW+Wf66LmX9NByhJT5O0\nunI7Ir6S3L9a0mER8c/taGAR+cQ3yz8Xk82/ucxB/SNwbsr9v00ee3lLWtQDarOFXPPLLJ+qkyYs\nf+YSoA6MiM21d0bEZknLW9aiHuET38xsYeYyB1VvW3dPrpiZWUvNpQe1QdJbIuKT1XdK+m/AxtY2\nq7hczt+s+Px3nA9zCVDvBK6S9Dp+F5BGgb2AV7a6YUXkcv5mxee/4/xoeogvIu6NiOcC7wduTb7e\nHxEnRMQ97WlesdQr529mxeC/4/yYy466+wBvBZ4GbAY+FRF7FpvrY16ga1Z8WX+vpalpVoyt9ZBf\nB81liO9zwAzwQ+AlwDMoD/v1teqx6kUSO6N2W0Iv0DUrkoOHhyhlBKnAQ36dNJcsvqMi4k8j4hPA\nq4E/bFObCqMyVl2amiYgNTh5ga5ZsaQttK/lIb/OmEsPanfZ34jYIakNzSmWrKKwAxK7IjwUYFZA\ntQvt97zsLPOQX/vNJUA9W9Jvk+8FDCW3BUREPL7lrcu5rLHqXRHcMn5Kh1tjZq1Suy2Hh/y6Yy5Z\nfAMR8fjk63ERsbjq+74LTuCisGb9wEN+3eMt3xfgnJOPYHDR7KHOwUXynJNZD6ktKpulMuTnreNb\nx1u+L1TtGeupObOe4yG/7nAPagEuWLeVmZ2zp1Bndoa7+mY9zEN+neMe1DxNTJYyr6K8MNesdzWb\n5efPgYUrRA9K0vsklSRtSr5e2s32VNY/ZXGShFlvW7NqhGvHTuKW8VMYcbJU2xQiQCUujIiVydfX\nu9mQrPVP4IW5Zv0mawftE49cyurx9U6cWAAP8c1Dva67t4w26y9pO2ifeORSrtxYckX0BSpSgHq7\npDcAG4B3RcS2tIMknQWcBbBs2bK2NOQJQ4NMTc/scf/I8JBPPrM+VLuD9urx9ZkV0f0Z0bzcDPFJ\nukbSDSlfpwL/AhwGrATuBv5P1vNExEURMRoRo0uXLm15OycmSzz82J5F3L3+ycwqvLNBa+SmBxUR\nL2zmOEmfBL7W5uZkSkstB9hvn8W+MjIzILsi+iLJ9fvmIDc9qHokHVR185XADd1qS9YV0NT2PYf8\nzKw/Za2V2hkxazGvEyfqK0SAAj4sabOknwEnAmd3qyGuv2dmjdSWRxpI2f3Bi3kby80QXz0R8fpu\nt6HinJOP4Nwvb541AerUcjOrVZ04sWJsbeoxnpOqryg9qNyoXBkNDw3uvm+fQb+NZpbNIy/z40/W\neXp0x67d32/bPuPxZDPLlLWY1yMv9TlAzUNaJQmPJ5tZlto5qZHhIS/qb0Ih5qDyxmsczGyuahfz\nTkyWWD2+fnf1Caed78k9qHnweLKZLUSl4HQpqYbutPN0DlDz4PFkM1sITxM0x0N885BWHNLdczNr\nlqcJmuMANU+148lmZs3KKoXkaYLZPMQ3T5UJTu/1YmZz5WmC5rgHNQ+VCU7v9WJm85G1h9QF67Zy\n9mWbPG2QcICah3oTnP1+QplZc6qnCXzRm85DfPPgCU4zayVn9aVzgJoHr4Mys1byRW86B6h58ASn\nmbVS1sXt8JLBvk7GcoCahzWrRnjVsSO793gZkHjVsU47N7P5SbvoHRwQDz2yo6+rTThAzcPEZIkr\nN5bYGeWt33dGcOXGUl+dOGbWOmnFZPfdazEzu2LWcf02L+UsvnlwFp+ZtVrt4n9vcuge1Lx4QtPM\n2s3JWA5Qc1KpHhEZj/fTiWNm7eVkLA/xNa12IV2tfjtxzKy9XJTaAappafNOFSN9eOKYWfv1e1Fq\nB6gmZc0vCbh27KTONsbMrA84QDXJ5fHNLA8mJkt9M+znJIkmecLSzLotbav4sy/bxPIerTThHlST\nPGFpZt2WNhdeySruxQroDlBz0O8TlmbWXY3WWvZawQAHqCb005ivmeVX1lx4tV4qGOA5qAYmJkuc\n86XrZ435nvOl63turNfM8i9tLrxWLyVu5SpASfpjSVsk7ZI0WvPYuZJulLRV0smdatP7rt6yR8HG\nmV3B+67e0qkmmJkBs4vKQnmZS7VeS9zK2xDfDcBpwCeq75R0FHAGcDRwMHCNpKdHRPrK2Raamp6Z\n0/1mZu1Uu1V8L08/5CpARcQvAKTa6wJOBS6NiEeBWyTdCBwH/KizLZxtYrLUUyeDmRVLbeJWpV5o\nrwSsXA3x1TEC3FF1+87kvrbbf8lg5mP9tC+LmeVb2hqpom9w2PEAJekaSTekfJ3aouc/S9IGSRvu\nv//+BT/fKc86KPOxXsqWMbNiq7dPXVF1fIgvIl44jx8rAYdW3T4kuS/t+S8CLgIYHR3N2hmjKe+Z\n2Mwl192e+XgvZcuYWbH14j51RRniuxo4Q9LeklYAhwP/0c4XnJgsccl1t2fu/dRr2TJmVmy9uMFh\nrgKUpFdKuhM4AVgraR1ARGwBLgd+DnwTeFs7M/gmJku86/LrM4MTwPmnHVPoyUcz6y29WC80b1l8\nVwFXZTx2HnBeu9tQmWjcGdnhaWR4yMHJzHKlF+uF5ipA5UG9jQmhvDCuyFckZta70tZInX3ZpsIG\nKweoGvUmFAW87vhlhftPNrP+UhkJqlxsF7XSea7moPIga0JxQOLC01fyoTXHdLhFZmZz0ysp5w5Q\nNbKKMT5+yJ1NMyuGXkk5d4CqUSnGODw0u4LEtu0zrmJuZoXQKynnDlApssZoXcXczIqgV1LOPW6V\nwVXMzayoeiXl3AHKzKwH1VY6LyIHqAyLBLsWVMnPzCw/irh3lOegMtQLTu+Z2Ny5hpiZLVBRt+Jw\ngEoxMVnaYyvlal/88R11HjUzy5eirotygEpxwbqtdQvF1qvTZ2aWN1nrn0pT06wYW8vq8fW57E05\nQKVotJhtYM8t6c3Mcqve+qc8D/k5QKVotJjttf/l0LqPm5nlSVaFnGrTMzt552WbctWbcoBK0Wgx\nm+vxmVmRVCrkjAwP1Z1fh3z1phygUqxZNcL+SwZTHxspWKkQMzMof65dO3YSt4yf0vBzLC8JFA5Q\nGd778qN7olSImVmtZob88lBY1gt1M/RKqRAzs1rVn2+ljEA0vGSQ1ePru/r5p+jhlOnR0dHYsGFD\nt5thZpZbtZsbAgwOCKJcILtiaHCA8087piVBStLGiBhtdJyH+MzM+lhtAsXI8BD77rV4VnCC7sxL\neYjPzKzP1RaWXTG2NvW4Ts9LuQdlZmaz5GXDQwcoMzObJS8bHnqIz8zMZslLFrMDlJmZ7SEPGx46\nQJmZWVM6vemhA5SZmTVUu16qUrMPaFuQcpKEmZk11I1NDx2gzMysoaw1UO1cG5WrACXpjyVtkbRL\n0mjV/cslTUvalHz9azfbaWbWb7qxNipXAQq4ATgN+EHKYzdFxMrk660dbpeZWV/rxtqoXCVJRMQv\nAOQt1c3McqUba6NyFaAaWCFpEvgt8J6I+GHaQZLOAs4CWLZsWQebZ2bW2zq9NqrjAUrSNcDvpTz0\n1xHxlYwfuxtYFhG/lnQsMCHp6Ij4be2BEXERcBGUt9toVbvNzKyzOh6gIuKF8/iZR4FHk+83SroJ\neDrgzZ7MzHpU3pIkUklaKmkg+f6pwOHAzd1tlZmZtVOuApSkV0q6EzgBWCtpXfLQHwI/k7QJuAJ4\na0Q80K12mplZ++UqSSIirgKuSrn/SuDKzrfIzMy6RRG9m0cg6X7gtgU8xZOAX7WoOUXm96HM74Pf\ngwq/D2XzfR+eEhFLGx3U0wFqoSRtiIjRxkf2Nr8PZX4f/B5U+H0oa/f7kKs5KDMzswoHKDMzyyUH\nqPou6nYDcsLvQ5nfB78HFX4fytr6PngOyszMcsk9KDMzyyUHKDMzyyUHqDokvUtSSHpScluS/knS\njZJ+Juk53W5ju0i6QNIvk9/zKknDVY+dm7wHWyWd3M12doKkFye/642Sxrrdnk6RdKik70r6ebKR\n6J8n9x8g6duS/jP5d/9ut7UTJA1ImpT0teT2Ckk/Ts6LyyTt1e02tpukYUlXJJ8Nv5B0QjvPBweo\nDJIOBV4E3F5190so1wE8nPKWHv/ShaZ1yreBZ0bEs4D/B5wLIOko4AzgaODFwP+t1EnsRcnv9nHK\n//dHAa9N3oN+sAN4V0QcBRwPvC353ceA70TE4cB3ktv94M+BX1Td/nvgwoh4GrANeHNXWtVZHwW+\nGRFHAs+m/H607XxwgMp2IfBXQHUWyanA56PsOmBY0kFdaV2bRcS3ImJHcvM64JDk+1OBSyPi0Yi4\nBbgROK4bbeyQ44AbI+LmiHgMuJTye9DzIuLuiPhp8v2DlD+MRij//p9LDvscsKY7LewcSYcApwAX\nJ7cFnES5Nij0wfsg6QmU66J+CiAiHouIKdp4PjhApZB0KlCKiOtrHhoB7qi6fWdyX697E/CN5Pt+\new/67fdNJWk5sAr4MXBgRNydPHQPcGCXmtVJ/0j5gnVXcvuJwFTVRVw/nBcrgPuBzyRDnRdL2pc2\nng+5KhbbSfU2TgT+N+XhvZ7WzOaRkv6a8lDPJZ1sm+WHpP0oF2t+Z0T8ttx5KIuIkNTTa1UkvQy4\nL9mL7vndbk8XLQaeA7wjIn4s6aPUDOe1+nzo2wCVtXGipGMoXylcn/whHgL8VNJxQAk4tOrwQ5L7\nCqnR5pGS/gx4GfCC+N2CuZ56D5rQb7/vLJIGKQenSyLiy8nd90o6KCLuToa47+teCztiNfAKSS8F\n9gEeT3l8mfgZAAADTUlEQVQuZljS4qQX1Q/nxZ3AnRHx4+T2FZQDVNvOBw/x1YiIzRHx5IhYHhHL\nKf+nPCci7gGuBt6QZPMdD/ymqmvbUyS9mPKQxisiYnvVQ1cDZ0jaW9IKygkj/9GNNnbIT4DDk4yt\nvSgniFzd5TZ1RDLP8ingFxHxkaqHrgbOTL4/E/hKp9vWSRFxbkQcknwenAGsj4jXAd8FXp0c1g/v\nwz3AHZKOSO56AfBz2ng+9G0Pap6+DryUcmLAduCN3W1OW30M2Bv4dtKTvC4i3hoRWyRdTvnE3AG8\nLSJ2drGdbRUROyS9HVgHDACfjogtXW5Wp6wGXg9sTjYLhfLw9zhwuaQ3U97O5jVdal+3vRu4VNKH\ngEmS5IEe9w7gkuRi7WbKn4GLaNP54FJHZmaWSx7iMzOzXHKAMjOzXHKAMjOzXHKAMjOzXHKAMjOz\nXHKAMjOzXHKAMjOzXHKAMushkv5Z0k8l/X6322K2UA5QZj0iqSz9ZOC/U66haFZoDlBmHSBpp6RN\nkm6Q9CVJS5L7f0/SpZJukrRR0tclPb2J51suabqqBBER8TBwEPA94J+S44aS131Myc7QZkXhAGXW\nGdMRsTIingk8Brw1KcZ6FfC9iDgsIo6lvHNxs/vp3BQRKys3JD0RWAI8SLlOIhExnRxzVwt/F7OO\ncIAy67wfAk8DTgRmIuJfKw9ExPUR8cN5Pu97gH8AtgBHL7iVZl3mAGXWQZIWAy8BNgPPBDa26HmX\nA88FLqO8NbsDlBWeA5RZZwwl80UbgNtpsDWDpKdK+pSkK5p8/g8BH0g2lnSAsp7g/aDMOmO6er4I\nQNIWfrfh3SwRcTPw5mYClKSVwGnA8yR9nPKur5sX3mSz7nIPyqx71gN7SzqrcoekZ0n6gzk+z99T\n3vm4sgv0s3EPynqAA5RZlyTDca8EXpikmW8BzgfuafY5JJ0ELImIa6qe915gP0kHtLrNZp3kHXXN\ncihJGT8P+CPg4og4v+bx5cDXkrT1Zp7vVmA0In7V2paatY8DlFkBSToU+Hfg17VzWzXHDQE/ApYC\nx0TEAx1qotmCOUCZmVkueQ7KzMxyyQHKzMxyyQHKzMxyyQHKzMxyyQHKzMxyyQHKzMxyyQHKzMxy\n6f8DGCoGp9BqInsAAAAASUVORK5CYII=\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x7f1bc8372290>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"f, ax = plt.subplots()\n",
"ax.plot(reduced_data[:, 0], reduced_data[:, 1], 'o')\n",
"ax.set(xlabel=r'PC$_1$ [$\\AA$]', ylabel=r'PC$_2$ [$\\AA$]', title='PCA of ADK')\n",
"f.tight_layout()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Convenience functions to create a new analysis"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"import numpy as np\n",
"import MDAnalysis as mda\n",
"from MDAnalysisTests.datafiles import PSF, DCD, PDB, XTC\n",
"\n",
"u = mda.Universe(PSF, DCD)\n",
"backbone = u.select_atoms('backbone and protein')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We added a new analysis class `AnalysisFromFunction` to make it easier to write a new analysis. Now code like this with a handwritten loop."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 0.11819594, -0.00500871, -0.32482931],\n",
" [ 0.11087458, 0.00267215, -0.29074168],\n",
" [ 0.09788045, 0.02156858, -0.29609087],\n",
" [ 0.02365893, 0.04328015, -0.27346295],\n",
" [-0.04413487, -0.03304375, -0.19893967]], dtype=float32)"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"result = []\n",
"for ts in u.trajectory[:10:2]:\n",
" result.append(backbone.centroid())\n",
"result = np.asarray(result)\n",
"result"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Can now be converted into this."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 0.11819594, -0.00500871, -0.32482931],\n",
" [ 0.11087458, 0.00267215, -0.29074168],\n",
" [ 0.09788045, 0.02156858, -0.29609087],\n",
" [ 0.02365893, 0.04328015, -0.27346295],\n",
" [-0.04413487, -0.03304375, -0.19893967]], dtype=float32)"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"import MDAnalysis.analysis.base\n",
"\n",
"cog = mda.analysis.base.AnalysisFromFunction(\n",
" lambda ag: ag.centroid(), backbone, stop=10, step=2).run()\n",
"cog.results"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"If you have a specific observable that you want to calculate several times you can also create a new analysis class with `analysis_class` like this."
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Step 49/49 [100.0%]\n"
]
}
],
"source": [
"def cog(ag):\n",
" return ag.center_of_geometry()\n",
"\n",
"COG = mda.analysis.base.analysis_class(cog)\n",
"\n",
"cog_results = COG(backbone, step=2, verbose=True).run()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Calculating the centroid for another simulation now becomes trivial."
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[ 60.18655014, 52.12106323, 28.70299911],\n",
" [ 60.20291519, 51.26913071, 28.00538635],\n",
" [ 60.38302612, 50.37888718, 28.21063995],\n",
" [ 60.63026047, 50.09779358, 27.16624451],\n",
" [ 59.51187134, 50.36880875, 26.06419945],\n",
" [ 58.09186935, 53.2186203 , 26.51293564],\n",
" [ 57.44642639, 51.7155838 , 26.66236305],\n",
" [ 57.62544632, 51.99805069, 25.02329445],\n",
" [ 56.43128204, 52.2718811 , 27.21855164],\n",
" [ 56.57172775, 52.93976974, 28.06817627]], dtype=float32)"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"u = mda.Universe(PDB, XTC)\n",
"backbone = u.select_atoms('backbone and protein')\n",
"\n",
"COG(backbone).run().results"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Using this classes you will automatically benefit from any performance improvements in MDAnalysis."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Memory Reader: Reading trajectories from memory\n",
"\n",
"MDAnalysis typically reads trajectories from files on-demand, so that it can efficiently deal with large trajectories - even those that do not fit in memory. However, in some cases, both for convenience and for efficiency, it can be an advantage to work with trajectories directly in memory. In this release, we have introduced a MemoryReader, which makes this possible. This Reader has been originally implemented in the encore package.\n",
"\n",
"The MemoryReader works with numpy arrays, using the same format as that used by for instance DCDReader.timeseries(). You can create a Universe directly from such an array. Here we create a random trajectory of 100 frames for the PSF topology."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"import numpy as np\n",
"from MDAnalysisTests.datafiles import PSF\n",
"from MDAnalysis.coordinates.memory import MemoryReader"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [],
"source": [
"universe = mda.Universe(PSF)\n",
"coordinates = np.random.uniform(size=(100, universe.atoms.n_atoms,\n",
" 3)).cumsum(0)\n",
"\n",
"universe2 = mda.Universe(PSF, coordinates, format=MemoryReader)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Certain operations can be speeded up by moving a trajectory to memory, and we have therefore added functionality to directly transfer any existing trajectory to a MemoryReader using Universe.transfer_to_memory:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [],
"source": [
"import MDAnalysis as mda\n",
"from MDAnalysisTests.datafiles import PDB, XTC\n",
"\n",
"universe = mda.Universe(PDB, XTC)"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10 loops, best of 3: 96.1 ms per loop\n"
]
}
],
"source": [
"%%timeit\n",
"bb = universe.select_atoms('protein and backbone')\n",
"cog_positions = np.empty((universe.trajectory.n_frames, 3))\n",
"for i, ts in enumerate(universe.trajectory):\n",
" cog_positions[i] = bb.centroid()"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"universe.transfer_to_memory()"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10 loops, best of 3: 53.3 ms per loop\n"
]
}
],
"source": [
"%%timeit\n",
"bb = universe.select_atoms('protein and backbone')\n",
"cog_positions = np.empty((universe.trajectory.n_frames, 3))\n",
"for i, ts in enumerate(universe.trajectory):\n",
" cog_positions[i] = bb.centroid()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"You can also directly copy a trajectory to memory during the universe initialization."
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"universe = mda.Universe(PDB, XTC, in_memory=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Likewise, the AlignTraj class in the analysis/align.py module also has an in_memory flag, allowing it to do in-place alignments."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2 (mda2)",
"language": "python",
"name": "mda2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.12"
},
"toc": {
"colors": {
"hover_highlight": "#DAA520",
"running_highlight": "#FF0000",
"selected_highlight": "#FFD700"
},
"moveMenuLeft": true,
"nav_menu": {
"height": "99px",
"width": "252px"
},
"navigate_menu": true,
"number_sections": true,
"sideBar": true,
"threshold": 4,
"toc_cell": false,
"toc_section_display": "block",
"toc_window_display": false,
"widenNotebook": false
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment