{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Nilearn statistical maps. Plotly interactive brain slice selection and plotting ##" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Recently a new `nilearn` [https://github.com/nilearn/nilearn](https://github.com/nilearn/nilearn) version was released. \n", "The interactive visualization of statistical map slices, as a new feature introduced in this version, suggested a few Plotly experiments intended to visualize more info on slices during cut point selection, brain slice plotting or 3d animation." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`nilearn` works with NIfTI images, a common format in neuroimaging (brain imaging). The NIfTI images can have at least 3 dimesions.\n", "The first 3 dimensions are spatial, the fourh one, if exists, is a temporal one, and other dimensions are user defined.\n", "\n", "A 3d image array has as basic elements the voxels. A voxel index, (i, j, k), points out its position in the image.\n", "But the voxel coordinate system does not give any information on the physical dimensions (units) or\n", "on how the head (brain) relates to the voxel indices. That's why the volumetric description of the imaging data as a 3d array is complemented with a $4\\times 4$ affine transformation, that maps a voxel index, (i, j, k; 1), to MNI (Montreal Neurological Institute) spatial coordinates, $(x,y,z;1)$.\n", "\n", "The MNI coordinate system is a right coordinate system (i.e the cross product of the unit vectors of Ox, Oy axes is the unit vector for the Oz axis), that has the origin at an anatomical point called, *anterior commissure* [http://imaging.mrc-cbu.cam.ac.uk/imaging/FindingCommissures](http://imaging.mrc-cbu.cam.ac.uk/imaging/FindingCommissures) and the ] axes are oriented like in the image below: " ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAMgAAADcCAYAAAA1H+4TAAAABGdBTUEAALGPC/xhBQAAACBjSFJNAAB6JgAAgIQAAPoAAACA6AAAdTAAAOpgAAA6mAAAF3CculE8AAAABmJLR0QA/wD/AP+gvaeTAAAAB3RJTUUH4wQEAC8VWScunQAAS6FJREFUeNrtfXdcVUf6/nNuo/cqiCBiwV5Bxa4Ea2LsJbHFkrY/NZuYvrvJN9kkm5iyycZgicYSo8aWxBZbsLdYQIOoCAoizSsd7j3l/f0BczwXsaDIveB58iHgLefMmZln5m3zvoAKFSruCI6IlgGwAyAB4KzdIBUqbAkcEZkB6K3dEBUqbBE6ACYAWqg7iAoVt0FX8aOp+FGhQoUCKilUqLgLVIKoUHEXqARRoeIuUAmiQsVdoBLERkBEICJrN0NFJagEsQFIkiSTQyWJbUEliA1Ao9GgrKwMBQUF4DjVFWVLUAliZUiSBAA4duwYYmNjLV5TYX2oBLEymEh14MABfP/99ygsLIRGo4Eoiqq4ZQNQCWJFEBG0Wi2ICIcOHUJSUhL2798PABBF0drNUwGVIFYFE6USExPx559/AgDWr18PAKouYiNQCWJFMILExcUhJycHBoMBO3fuxJUrV6DX61URywagEsRKYOKVyWTC999/L/tB0tLS8PPPP8ufUWFd6KzdgMcVRASNRgOj0YiWLVsiKCgIOTk5AACj0QhJksBBAlUYtCRJAAcOEjgIEKHVaKHR6MFBXeUeJTgiKgVgb+2GPG5QOgc5joNGo8Gvv/6Knj17wtnZufw1EAARnEYDcJZn2iQQBCLoOI1KkEcIdQexEpRKONNFBEGARqOBTqeDJBHMggg7PQdj9jXsP3IGTi6ucNSKyMgrRkiTMLRr3QoEAoFTT7o9IqiLj5XAdg0mamk0GnAcp9hVAI0GAGfAiWPHsXfPHnh7e+PfH/4bo4cPx7lzF6EHQJJQTg5VXXkkUHcQK4P5QYBbpCn/G9BpCICE5q07oM8Tw5Bw5hR2/3EQQ56ajrGjh8NMArQyO9Q95FFA3UFsDXLQIsARD1GUENg4DDqY8f5br4Dn7PH3V1+BvQaQJALHQd09HiFUgtgaKnQTjgOI0wOSGToA69etwy+7j2DKy3PQt0c4srNzoIMGGmjLP6xuII8EKkFsFhwk0kCr1yMj+SzefONNeDdsgnmvv4Ijx49h3do10Gm1IFKZ8SihEsRGcctJqMefpxPg7ReAkcOH4uSRo/huwUK4u7tX+pyKRwHVD2IDICJwHIfNmzejX79+cHFxUbwpQhRFaHV68IKIgpJSONvbwc5gsHazHwuoVqxHgMqnA++2ynMcB0mSoNFoIEmSHMUriuUOQg4aaHVaAIBep4OXqwtU1B5UgjwkiMjigJOlqfb+9AP2eUdHR9jZ2QEoN/8yiEQgSQJHEgAJnEYHjtNCDfh99FAJcp9gJFASgjn4lJMZAMxmM/Lz83Hjxg1kZGQgNzcXhYWFKCkpgdFoRHFxMcxms+w15zgOrq6uSE9Px9WrV+Hr6ws/Pz8EBATA3d0dzi4uFffQWrSH5wXZwQjgtnaoeHioOsg9wKJsmcdbuSsQEQoKCnD+/HkkJSUhISEBV65cwZUrV3Djxg3Y2dnJBLK3t4eHhwfs7Ozg4OAAURQhSRK0Wi0kSYLZbEZGRgYEQQARQRRF8DwPrVaLoKAgNGrUCOHh4Wjdug3Cw1sgICBA3nkqk1Y9S1JzUAmigDL1DlOc2SRkr12+fBknT57EsWPHcPbsWaSlpcFkMsHb2xshISEIDg5GmzZtEBAQAF9fX/j4+MDNzQ0Gg+GeE/e3335Dr169YGdnh+zsbGRlZSE7OxuJiYlISkpCcnIyrl+/DkEQEBQUhM6dOyMqKgrdunWDj4+PfB0WCFm57RzHqeSpJlSCVIApx2xVZ5MrPz8fBw4cwM6dO3HixAlcu3YNDg4OaNasGTp16oT27dujRYsWaNy4MXS6e0usjISVdyKNRlO1FavSd2/cuIGzZ8/i6NGjMkk5jkObNm0QHR2NmJgYBAcHy98RBEH+m4mEKu4fKkFwixxMhi8pKUFcXBzWrl2LI0eOQBRFtGzZEr1790ZkZCRatWoFNzc3i2tUpaxX9bsqMCvW5s2b0bdvX7i6ulpci12b6SvK141GI06cOIFt27YhLi4OBQUFaNeuHUaPHo2BAwfCw8OjymdUcX947AiiXMElSbIQoy5duoTly5fjl19+QV5eHiIiIjBq1Cj07dvXQoQByie1UpSpLL5UR5S5qx9EcT8GURTl+yknfGlpKQ4fPoyff/4Zu3fvBhFh0KBBmDZtGtq1aydfp7L4pYpdd8ZjRxCmBCutT0ePHsX//vc/7NmzByEhIRgzZgyefvppBAUFyd9T6ic1LcvfD0Hu9t2qjAi5ubnYtm0bfvzxR5w5cwbdu3fHSy+9hF69ekGr1cr9wEimkuQOIKJSekwgSRLxPE+SJBER0YkTJ2jkyJHk7e1Nw4cPp507d5LZbJY/bzabied5EkVR/s6jahcR0aZNm6igoOCBryMIAgmCYPGMoijSkSNHaOLEieTv708xMTG0Y8cO+X32fCqqxmOhsRGRvGLqdDqkpKRg6tSpGDRoEHiex+bNm7Fx40YMGDAAWq0WolgR3qHVyqtrXVhhmXGBnTFhekdkZCRWrlyJPXv2oHHjxpgyZQqGDRuGY8eOyXqNIAhqRscqUO8JwiaKTqcDz/OYP38+evbsiatXr2LNmjXYtGkTunfvLivCbIKxyVZXyMHA2sueg4WwEBFatGiBBQsWYOfOnfDy8sJTTz2F6dOnIzU1FTqd7jZDgwrUbxFLkiRZfDh58iRFRUVRUFAQLV26VBalzGazhVhlrXYSPbyIdS8w8YvhwIED1LdvXwoICKCvv/5a7odHLVLWJdQrgjBCiKJIZrNZ/vuLL74gT09PmjBhAl27do2IbsnrkiRZfTLUFkFY/7BnJyrXQRYtWkQNGzak/v37U0JCgtw/qn5SD3UQqtA39Ho9jEYjRo8ejU8//RT/+9//sGrVKjRo0EA2k9Yl/aImUFn0EkURGo0G06dPx+HDh9GgQQP069cP8+fPB1DuWGS62+OKekkQg8GA+Ph49OrVCzdu3MDBgwcxbtw48Dwv6xmPCynuBuYLEQQBgYGBWLFiBb755ht8+eWXePrpp5GRkQGDwfBYZ5qv8wShCsWSWZ50Oh3Wrl2Lvn37omfPntixYwdCQkLk9x63XeNOYH3AdlKqMGaMGTMG+/btg8lkQq9evfD7779X5OmSHkui1HmCcBwnmzN1Oh2+/vprzJw5E6+++ipiY2Oh0+lkk+3jToqqoLR6aTQa8DyPxo0bY8uWLZgyZQrGjRuHTz75RO6/x87KRXVcSRdFUbbMvPfee+Tq6kpr164lonIFVOk0s1XUlpJ+P+2QJMlCif/tt9+oQYMGNGXKFCoqKiIikt97HFBndxAWC8WC+ObPn4/PPvsMq1atwujRoyEIgoWjT8W9oRS52G4yZMgQ7N69G2fOnMHgwYORlpYGrVYLnucfC3GrzhIEgKxXLFmyBP/+97+xYsUKDB06FDzP10knny2B4zjZuRoeHo4dO3bA1dUVPXv2xNGjR6HX6y1C6esr6ixBRFGEXq/Hli1bMG/ePHz55Zd46qmn5J1DJcbDge0OjCQ+Pj5Yt24dBg0ahGHDhmH79u3Q6/X1XnGvkwRh5IiPj8f06dMxb948PPvss7IyrppxHx5KKxcjgsFgwIIFCzBnzhyMHj0aa9eulYsA1Vflvc4lbaCKykw3b97E5MmTMWzYMLz++uvy+W2VGDUPqgilpwpT8FtvvQUHBwc899xzkCRJ9jHVx9OKdYogpDjsNHfuXOh0Onz++efy+yo5Hi1YeQae5zF37lw4Ojpi5syZ4DgOY8eOleub1Cei1CmCMKV82bJl2Lp1K7Zv3w5nZ2dZ5FLxaKBceJRhKrNmzUJRURGmT58ONzc3DBw4EGazGXq9vt4sVnWGIMyce/HiRbz77rt4++230bFjR/A8r5KjlsEci4Ig4O9//zt4nse4cePkEnJMF6wPsHmCkOKYKwC89tpraNOmDWbPni3vKCpqH0yBFwQBb7zxBm7cuIGxY8diz549aNGihRwIWdd3kjoxu1h07rp163D06FHs2rULQM2fDVdRPTDFXZIkfPrpp8jNzcWoUaOwe/du+Pn5gef52zKx1DXYvDZFRNDr9SgqKsKHH36IadOmoVWrVvIKpcK6UE7+hQsXIiAgAM888wzMZnO92N1tfoYx+/qiRYtQUlKCOXPmlDdcJYfVwcjBwlL0ej1++OEHXL16Fa+88ooc3KgseV3XYNOzjCqSLOTk5GDhwoX429/+Bh8fHwiCUKe37foENg46nQ6CIKBBgwZYtWoV1q1bhyVLlshxW3UVNk8QAFi8eDEcHR0xdepUOdmbCtsDI0Pnzp3xxRdfYN68eThy5Ajs7OzqrKfdJglCigwjRqMRK1aswMyZM+Hs7Cy/rsI2wWK3JkyYgGeeeQYzZszAzZs35Sz2pEjAVxdgkzONFDmdNm7cCAAYPXo0AKjhJDaKqk4ofvTRR3BycpL1xrpEDAabJAgLkBMEAT/88AOefvppeHp6qrpHHQHzjzg6OmLx4sXYtm0bli9fLiflq0uwSYKwTjxw4ABSU1PxzDPPlDdWFa3qBNhZEkEQ0Lp1a7z77rt466235AR1dekciU3OOKbQrVq1Cp07d0Z4eLicqkdF3QATtURRxN/+9je0adMGc+fOld+rK+KWzRGEKtL2sHQ9I0eOlF9XCVK3oByvr776CgcOHMDSpUvlLCl1ATZJEAD4448/QESIiYkBoIay11WwUgvNmjXDm2++iX/961/IzMy0sGrZMmyOIAw7d+5E+/bt4e3trZp26zhYZO+cOXPQqFEjvPvuuwBQJ3YRm5p17ORacXExjh8/jgEDBgCoGx2p4s5gucs0Gg0++OADrF+/HnFxcXLOMluGzRBE6fs4d+4c8vPzERUVVd7IehA2/biDedl79+6NESNG4J///GedCDi1mdYpCXDo0CEEBASgadOmt72nou6CkeGf//wnLl26hB9//FE+nWirsBmCKHHy5Em0b99eNhOqBKn7YOIzq/H+wgsv4LPPPkNJSYmssNsibIYg7EhtcXExEhMT0bFjR2s3SUUNonJV3pdeegllZWVYtGgRANhsfi2bIQhDeno6CgoK0KpVK2s3RcUjAFPY3d3d8fLLL+Obb76B0WiEXq+3yV3EZgjCVo+kpCTo9XqEhoYCUPWP+gimi0ydOhV2dnb4/vvvAcAmdRGbI0hKSgr8/Pzg5eWlnv2oxxAEAc7Ozpg5cyaWLFmC/Px8GAwGmxOzbIYgDOfOnUNISAgANbykvkI5plOmTIEoivjuu+8A2J7Py2YIwrbdnJwc+Pr6ArC9zlJRM2AWLZ7n4erqikmTJuGbb75BQUGBzYXE2wRBlCbA/Px8NGrUCICqf9RXsMRzLOvJtGnTQET46aefANiWRcsmCMJQVFQEo9GI4OBgACpB6js4jgPP8wgICMDkyZNlk68tZWW0OYKYTCZ4enpauykqaglMtJ45cybS0tKwe/dum3Ic2gRB2HZaWloKURTh6Oho7SapqCUwXSQ4OBi9e/fG4sWLrd0ky/ZZuwFKFBcXw9HREc7OztZuiopaglKMnjZtGg4ePIisrCw58YO1YVMEKS0tlRM2AKoO8riAKev9+vWDs7OznMnGFsQsmyJIWVmZXGMQUAnyuICFn+j1ejz11FNYt24dANtI0mH9Figbo9HA2dkZ9vb21m6KiloGWwyHDx+Os2fP4sKFC3JuX2vC5ghSWlqK0tJSazdFRS2DESQiIgJBQUH49ddfAVjfJ2JTBHFycoIkSTCZTADqZiY+FQ8GlmyO4zgMHDgQmzdvBmB9McumCKLT6epFTQkVDwa2iwwdOhSXL1/GlStXrG7NsimC2NnZobCwEEVFRdZuigorgFWs6tSpEzw8PPDHH38AsK4kYRMEYSuHk5MTiEgmSF3LBK7i4cHK7XXp0gW///67tZtjGwRhcHBwgJ2dHYxGo/yaaup9fKAc6wEDBuD06dMwmUxWFbNsiiAuLi7w9vZGVlYWAFVJfxzBfGBRUVEoKSlBYmIiAOvNBZsgCLN3GwwGuLu7Iy0tzdpNUvEIwAoj3W2yM6tVcHAwfHx8cOjQIfm71oBNEETZAS4uLsjNzbVqp6h4OCiJoPybnQPhOA5ms/mOTkB2YKpt27bYv3+/VZ/F5gjStGlTXL9+HQBUk28dhSRJsoOPkQIAcnNz8eeff2LlypXYsWOHRTbNqtCtWzckJSUBgNX0EJuZgUqCbNq0CaWlpXBwcFDPpdsg2Fgpx4XtFhqNBlqtVj7TkZGRgZSUFFy/fh0FBQXw8/NDZGQkAgMD7zi27LUOHTrAaDQiPT0dDRs2tMqz2gxB2CrTpEkTmEwmXLt2DWFhYXWCIGwlZInR6mM2ekYKRgLgVhgIqxvJdMn09HRcunQJmZmZcmnoiIgIBAYGyko4I9TdCNKkSRNotVqcPXsWDRs2hCiKtS5V2AxBGBo3bgyNRoOLFy8iLCzM2s25J9hAs0zlLFyiLhC7us8IlC9kV65cQXFxMVq0aAGtVouysjIkJycjKSlJ1jGcnZ3Rp08f+Pv7WywWyt3nTv3DiObm5gZfX18kJSVh4MCBj7eIpdFoIEkSfH19ERwcjNOnT2PQoEE2r6izDOXx8fEIDAyEl5dXvSo2qlSydTodvvzyS3z88ccwmUyYMGECBgwYgNzcXPj4+CAsLAyhoaEWB96UekZ1KhSzXbhdu3Y4ffo0AOv4xGxKDmBWjSZNmiA+Pt5qnVJdsGOjS5YsQWZmplyoUpIkq4drVwdsMrOdkCXuYxlIkpKS8OmnnyIrKwulpaWIjY1FWVkZZsyYgaeeegpt27aVa9mLoig/O9NLHgTBwcGy0cYai6XNEES55UZGRiIxMREFBQVyfI6tgllXOnXqhFGjRmHFihVISkqCTqeD2Wy2qRxPdwPbKdgz6XQ6cByH9PR0/PHHH1i1ahV2794NOzs7AIDJZIKbm5ucZFxJCDaOHLi7ilL3g7CwMBiNRvA8D51OV+tzwWZELOBWx3bq1AnFxcVISkpCly5dqpTnKyt51tpp2H1FUURoaCiee+45rFixAllZWejVqxfMZrPN6iNKy5OyKm1KSgouXLiA3NxcCIKAwMBA9O7dGw0bNkRISAjeeecd5Ofn4+9//zuaNWsGURRl8Yk9q1IZf5BnZ98JCgpCZmYm8vPz4e3trRIEAJo1awY/Pz8cO3ZMJghw+yqn7Hi2rbOBqW0rEkt85+HhgRdffBFLlizBzZs38dRTT1nUOLGmdYsRQkkKdg4jOTkZ58+fR35+PrRaLUJCQtClSxeLFEySJGHw4MHo27evnBWRWbGICIIgQKfTwZh9Hcf/PIngJi3RolljiGQGBz00D0AUDw8PSJJkQZDaXGxsRsQCblkvNBoNOnfujN27d5c3skKBlyRJtrHn5ubi3LlzOHv2rCyj6vV6+Xxzba80HMfJJki9Xo+ZM2eitLQUixcvhiAIFuSubVT2ZDPL07lz57BhwwYsXboUCQkJ8Pf3x7BhwzBx4kRERUXB09PTQkln4q6Dg4NMDuUOqtPpsG/fPvTr2wcDBw9Ft8huiF22GhrOAImXgGoMCbuum5sbnJycrBZdYVM7iLID+vXrh3nz5iE7Oxu+vr7yAG3atAnffvstkpOT4eDgAAAoKSmBj48PRo4cieeee07ODA/UvujFxAyTyYRx48Zh165diI2NxbPPPiuvho96F1GKn0pPdnFxMS5evIiLFy+ioKAAbm5uaN68OZ544gnZ8sR2AmVYiLIP2fMp+5ZZuIxGI86cicfXsd/iyoVTmDVrHj767GsMffJJBLo5VWv1Z59zd3eHnZ0dMjIyanUcGWyOIGwwu3fvDlEUsXv3bowfPx45OTl44YUXsG/fPkycOBH/+Mc/0LhxY3Ach8zMTPz2229YunQpvvnmG3z44Yd49tlnZcWRpRGqLXAcBwcHB0iShAEDBsDb2xsLFy7EhAkTEBQUBJPJBJ1O99ApNtnqzlZ2tsMqJ7XRaMSFCxeQlpYGo9GIoKAgtG/fHsHBwTAYDPK1lLvEvZxxlSc5+7dWq8WUKZPh4uKCnj06Y9PGbdh1Kge8IJXLKtVY/Nnz2NvbQ6fTITs7W369NmFzBGFilre3N1566SUEBASgoKAAgwcPhr29PY4cOXKbh71Bgwbo2LEj5s2bhwULFmD27NnYsGEDYmNj4evrK8vGtQ1WoJLVe1+zZg2io6PRtm1blJWV1QhBlA48Joqmp6fjwoULyMjIkPuyU6dOCA4OtrgnM+Mqd5mHaYubmytIAgS+DDodj8zcPAwbFI0QXxeIPA+t7v4XKiUJ7ezs5DwFtQ2bIwhwS+eYM2cOeJ7H0KFD4eHhgd9++w329vYwmUyyuREoH2ie56HX6zF37lwMHjwYkydPRt++fbFt2zY0atTIauEfWq0WJpMJDRs2xIwZM/DDDz8gOzsbAwYMqHaBIOXqySY1Mw5cv34d8fHxuHnzJgCgYcOG6N+/v4UnW5IkCIIgf++hScHaUvFbFCWQKEFnsMfmzZtgcPbEx//3dvm9OYK2mtIuGzN7e3urZbqxSYIwi5RGo8Gnn36KK1eu4MCBA7C3twfP89BqtYiLi8Ply5cRHByMbt26wcHBASaTCYIgoHnz5ti7dy/GjBmDmJgY7Ny5Ew0bNrTaTmJnZwdJkuDq6oqZM2fixx9/xLp16zB69Og7ioFKyx3TCbRarTypzWYzLl++jEuXLiEjIwMeHh5o0KABIiMj4e3tbXEd1pdKQ0JNgCr+00ADECrIoUfmtTScOnMBS5ctR6C/J8pKCqBzcAHhFpnuB0oTviAIFq/VFmySIEQEvV6P5ORkfPvtt/jvf/8Lb29v8DyP0tJSTJ8+HXv27EGDBg2Qn58PJycn/POf/8S4ceMgiiLMZjPs7e2xceNGPPXUU3j66aexd+9e2ctrjZ2E7YoajQZTp07Fhg0bsGjRIkyZMgV6vd6CEOw3+zzTFYqKinDlyhWcP38eRqMRnp6eaNSoEXr27Ak3Nzf5XkoDhdInUePjBA7lU16CRGboDDqUlRZi3cbNGDR4GNzdXLFj1x6IooTBMQPKn60aSrrSCMD6oNadxkRUSjYGQRCIiGjOnDkUHR1NREQmk4mIiF566SUKCwuj5ORkIiLKzc2lDz74gNzc3GjWrFlUVlZm8fn8/Hxq27YtTZo0iYiIJEmy2nNJkkRms5l4niciori4OPrvf/9L2dnZRES0fv16ys/Pt2hnfn4+HT16lH7++WeKjY2ln3/+mf766y8qKSmRryuKIvE8TzzPkyiKtfaMQsWPKJlJlIqouDiNpkx6kgAdObt7k1ZnIOjd6Je9h+R2Vuv6FfOgT58+9Mknn1i8VluwuR2EKhxPubm52LFjB959910AgMFgQGJiIjZv3oxly5YhNDQUkiTBw8MDb7/9NoYMGYIRI0ZgzJgxWL16NRwcHCAIAlxdXbF69Wr07t0bn3/+OV555RWrebaViblFUUSvXr3g6emJxYsXY9SoUXB1dYVGo4HRaMRff/2F9PR0FBUVwcvLCy1btkRYWJiFKKaMdbLGrsjhlh6i4fQ4f/4CCDqMnzIZnESQBDM8AoLRvm1rSKieeEQKA4Q1jw/YJEE4jsP+/fshSRL69+8vd9DWrVvRuHFj9OvXT/4cVfgc2rdvj7i4OAwYMACzZ8/GokWLZAW2ZcuW+M9//oNjx44BsG4AJCk82USE1q1bIyAgAN9++y1SU1ORk5NTYRFyQ9euXdGwYUNZb2DGCBYnZe0zJxXCIDhOB0ni0aF9Tyz7od9tnxNECbwgwVANLZ0qfCv5+fnIzMxEkyZNrPKMNkcQhn379qFDhw7w9fUFz/PQaDQ4efIkIiIiLDzubFUWBAFBQUHYvn07rl27BgAW8UVTpkzBlClTan33oEqmWGV4R1paGs6ePQuj0YiwsDCcPn0aN27cwMsvv2zxfWaYYHqMrYCrUNIBDTjOABIJIi+A4wSIGh0IHHTEQ6vRgNMZqnVtNr7x8fEoKSlBVFSU3Ie1CZsiCFtZAeD06dMYPnw4gFve2+TkZPTr10/+LAMTMYgIISEhaNy4sfw94FbEbW0RgxRnINhqD5SXd0hNTcWFCxdw/fp1eHh4oGHDhoiMjISvry/c3NyQl5eHFStWYMyYMTAYDBYWKHY9W4FGVtLLf3E6DhqdDoAOt8wCWsX/q4+4uDg0a9ZMXihr2+lrO72tgNFoxI0bN9CiRQsA5RO8pKREDga8ExiRqpJZa4McSscbm8j5+flITU3F5cuXkZOTAycnJ7Rp0wZ9+vSBq6srgFu7jMlkwogRI3DkyBEsWrQIzz77LNzc3OrVAaz7BRu/vXv3yruHNXZPmyIIW+XT0tJARBY7QUlJCTQajWzOvNtxzdqYTGxSkyIIkA1qQUEBzp8/j9TUVOTn58PT0xPNmjXDoEGDLGqfKI/nMkKVlpaid+/esvI+cuRIhISEyBHBtfV81gSTJHJzc3Hx4kW89dZbAKyjO9oUQRiMRiMEQbCw7TNZXPlva4AUjjuluMMsT6mpqZAkCZ6enmjbti2aNGlym+Wp8i5DiuA/5rNo06YN3N3dsX79evTr1w9t27YFz/PyZ+ozSdiic/r0adjb2yMiIgIAal28AmyUICaTCXq9vspKU9aYGErLEzMKsJQ2SUlJSE1NhVarhZ+fH7p164aQkBCLA0MPYo4VRRFBQUGYNm0aVqxYgZycHPTv318Ona/vBAGAXbt2ITQ0FC4uLlYz9dokQe4mazLl91FOEOWkVoZ9i6KIjIwMJCYm4urVqzAYDAgJCUFMTAwaNGhg4flVikQPMrBarRaCIMDFxQXPP/88Vq5cibVr12LMmDG3xVTVN7BnOnLkCHr27AnAehKDTRIEsOwQqjikw3Ec8vPzb3u/pu6n3CnYDmAymZCcnIzLly/j2rVrcHBwQPPmzfH000/Dy8vL4vtK8akmInXZGWyNRoPJkyfjl19+wZIlSzBx4kTY2dnBbDZbBG3WB7Dd0Wg04sqVK/jXv/4FwHq+K5skCBv80tJSWQ9xcnKCTqer8eI6yknNBqGwsBCXL19GSkoKjEYjnJycEB4ejt69e8PFxUX+Lgugq3wGoyagDNQjIpSUlODJJ5/E4cOHERsbi4kTJ8Lb27veJaljZ1pOnz4NIkKbNm0s+qO2YZME8fT0hEajQV5eHvz9/eVV1NPTEzdu3ABw7w5ju4FS/mdgZ9eVIkphYSESEhKQlpaGkpISuLq6omXLlggNDbVYpZVm5NoqV81xnBzJ3K1bN7i4uGDJkiUYOXIkwsLCLByJdRlKp+rRo0cRHBwsnw5VQ01wa6I1aNAA9vb2yMrKQosWLeQOatSoEVJTU+/rWspkZ8Ctic0mEhEhOzsbFy9exOXLl2XLU/v27REaGnrHmKdHGR17NzAyi6KI1q1bw9vbG2vXrkVkZCQiIyNly1pdVd6VzlUAOHjwoGy9UglSCV5eXjAYDEhOTkbv3r3lCern53dfWfaUVp5du3bh6tWriImJQWBgINLT03H+/HmkpaVBp9PB19cXPXr0QKNGjeSJz2KemCnWllZmjUYDk8kEf39/TJ06FStWrIDRaJSzUNaUhevW7kvQaB6tWZkRgC1mJpMJly9fxtSpUwFYN3bOpgjCYqz0ej1atGghZ1dkCA4Oxr59+2A2my3OBzCRSZleR6vVIjY2Fi+++CIkSULr1q3x3HPPwdXVFWFhYRg0aBD8/PxuSx0kd4wiRMSWwHEcDAYDJEmCs7MzZs2ahTVr1mD16tUYM2aM3IcPtssRJOIBgYNGp6+4BkEgAaKkhx1JIA0gcOULhhZUEW5SzbsoonTZLp+cnIzff/8dL7zwApKTk1FcXIy2bdtau7ttK+0PcKvzwsPDkZCQAOBWWa7mzZsjOztb1kNYKiD2HZYIwWw249y5c/j222/lJMhnz56Fg4MDpk2bhqioKPj7+8vXUAYTsh9bJAcDMx0zp+WECRNkz3tJSYmcYf5B+1+j1yDrWgo2r/sJ27dsRV5hKaABREkESAKH8onzoJOHOX3ZrvHzzz8jMjIS//d//4eioiIcP34czs7OCAsLq/ax5JqGzRGEdUb37t1x9epV2QkHlJ+z1mq1FiXatFot9Ho9iouLceLECaxevRpr1qzBxYsX0a5dOwDl8VAuLi5ymkxBEGQLlK2T4V59xfICx8TEoF27dli8eDGys7PlyN/qmMOJCFpOi3Nn4/Htgq8w/z8fYtDQoZg46UVkFpQAei0I5YGHWglANTmo1DMMBgOMRiNmzpyJ0aNHo3379ti1axecnZ1x5MgRtG3b1iJq21qwKRELsCye4ujoiLi4OISEhMBsNsPf3x+urq44f/48IiIikJeXh4sXLyIlJQWFhYXw8fFB69at0bRpU9jb2yMqKgpeXl44e/YsXnrpJXTu3BmSJMlWKaXvo3wyARqNdZTwB+0rZkwQBAFdu3aFn58ffvzxR0RHR6NVq1YQBMEimvluiwHHcTCVmmDmNXj9rQ/xj3f+iaFDBmP7to1Iu/Yegl2bQCSxWkIV61tl0j9RFLF69Wq8/vrrKCwsxMKFCzFlyhRZB2HlDpRjpJp5K8BxHHieh4uLC7p27YotW7Zg8uTJ0Ol0MBgM8PX1xc6dO+Hs7IzCwkL4+vqic+fOCA4OtjitZzab4enpiS+++EK+NvM+V8Ytc6/1sh8+LJjnvXHjxnj22WexfPlyZGRkIDo6Wl6177kSE0GvM6Bd+3bQcARIPBwcHTB20hS0bd4EkihAw3EQUX60/G7nn24tOiQTg+d5bN26FR999BEOHTqE8ePH49///jdCQkJko0hOTg4uX76Mzp07y21WRSxYrjRsIMePH4/4+HicOHEC+/btw08//QQAuHz5MiIiIjBx4kQMGjQIYWFh0Ol04HneohYFk9EFQYDZbL5tBSWUE9JsKkX8qT/x4Yef4HLqVbk9dQ3sIJaXlxdefPFFZGZmYu3atbIodq9M80QcUJFX4XrGebw46Rls2rIXnOiEsiITNJwGROX9JnIVEpZiF2ZlE5QRzjqdDnl5eViyZAkiIyMxbNgwCIKAuLg4/Pjjj2jUqJFFzquUlBSIooimTZvKz2TVsSArJm2QJIkEQZCTDTDwPE/Jycm0adMmCgwMpGeeeYYOHjxIBQUFtGHDBoqIiKDCwkIiKk/O8CBJCiRJIl6UKCs3h35eFUutAzzI4BJEJy+mEVH1Eww8bD8QEW3evJkKCgoe+jqiKMrt37ZtG3377bdyfwmCQJIk3bHPJMFMormESgpzac+W9TSkTzcCQP/6+kcSicgklJFEPElExAsCCYJQZV9lZmbSpk2baPr06eTl5UVarZaio6Np06ZNcmINnufl9rBkDAsWLKAuXbrc9kzWwiMTsUgR11TV65Vjlliep8TERBQWFoLjOHTs2BEzZszAzp07ERkZCa1Wi7Zt2yIjIwNnzpxBVFTUg4d+cxyIyk2lI0Y9jbWLFyH5VB4g6x+1v62z1Zf1U3WfSxmewgIaBw4ciKNHj2Lx4sUYP348/Pz8ZL2kilEDp9WC0+phr7VD38Ej0LpFKCK7R+PChaRyyxV3y7Cr4TSQpPITj5mZmTh9+jSOHTuGffv24dixY3KZtueeew4TJkyQjSYA7ni+JT4+Xt49lNn6rYUaJ4iSAIwclUnBOqSoqAiXLl3ChQsXUFRUBGdnZ4SGhqJZs2ZwcXEBx3Hw8PDA999/j7Vr12L8+PFo3LgxunXrhvfffx/bt2+XJ0PlOiHKiQZANotWnnR2BjtwnBkGOzuIgiCn6bCG1Ovt7V1j4Sus/wVBQGRkJLy8vLB27VrExMSgWbNmFsq7chEzFRlx+lwyQpu3gI+7G3yCQtGseRiiupTHRJFIgEaLchdieWdNnz4dP/74I3ieh5+fH7p06YKPPvoIvXr1Qnh4uBx0qTxnPm/ePEycOBHPPvushVh96tQpjBw50gq9fwdQDYtYLB/VjRs36JdffqHSUsvL37x5k44dO0bLly+npUuX0rZt2yg+Pl7edhlYriciojfeeIN69OghfyYtLY1CQkLo448/lu9ZOR8U+75S7GDbORGRRESCRCRKRET59ExMFNm5htLplMyKz1t3a69JMBEoOzubvvzySzp48CARlYs4ZrOZJEmSn/f8iT3k7migFu06U9yxk7Rt+w6aP38+FZSaSJSIJEEkkiQ5JxYR0ffff09ffvkl/fnnn7eJiCxfF1H52L/xxhtkMBgoPDyc4uLiLNp38+ZNat68OW3bto2Iaj8HVlWoMYIo5cjExETq2bMn+fr6UkpKCmVnZ1NcXBwtW7aMFi9eTNu2baOUlBQym823DaRSnmWTOT09nYKDg+m7776T39u8eTO5u7vTrl27iOgWSRhYWzZs2ECfffaZPFgWBCGi8m+UE8TgGkqnU65XtKX+EIT1LRFRcXExfffdd7R169bb+oRIIr70Jm3/bSO9+OJL9PH8L+noqQQqMfEkSgKJPF++qkjlfcdL4m2TWBRFMplMMvGIyon5ySefUIMGDcjV1ZU+++wzKi4ulseJte3UqVMUHBxMGRkZ5a2xsv5BVEMEUSp9GzZsoMDAQAJAer2eXn31VVq9ejXt3LmT0tLSLB66fOWyXNkrg60+CxYsIH9/fzmjIhHRO++8Q0FBQZSSkmLxWeXfS5YsobCwMDljIRsMiYh4SSJB4omkPBo3sBfZu4bQqUsZFZ+TbGKAWN8qf1ifMeOGUKEss9WaZW9kSjBbeEwmk9zX69ato+XLl8v9xD4jVKFw86JAomgiEkWSKggiEpFQsSgWFxdTSUnJbWRJSkqiefPmkbe3Nzk4ONDzzz9PV65cke/H2s/asHz5curYsaP8b1vo/4cmCNs5srKy6MUXXySUS/FkMBgIAL388styh7DPVyc9JhtkIqLo6GgaNWoUEd2yxjz99NPUrVu321Yk9p2LFy9SaGgonTp1Sn6fqJwgZkkgQSolIp5GDo4mJ7cAung1x+L6tQ3WP2zilJWVyStyTbWJ9cGOHTssLFxs1RcEUe5HSZJIoqrvqexnhitXrtDChQspOjqatFoteXh40MsvvywvbGy+KJ+DEeK1116jIUOGWLxmbTyUkk4VSm9hYSE++eQT7N27F61atUJpaSnKyspgNBpx+vRpFBUVwcnJySLZ2/1CmSzt008/RUxMDDZs2IARI0ZAkiQsXboU0dHRGD9+PNasWQM7Ozs5lxQRISAgAF5eXjh37hzat29v2X6JUFRYiFNHDyHx/EUUFxVg0y+/YvLYYfCsqAb1KK0odIejvcp6hpXzYEmShNzcXOTn56OgoABlZWVyfFpeXh6Ki4uRl5eHsrIymEwmFBcXo7S0FJIkoaysDL6+vvjvf/8Ld3d3PPHEEzh58qScPSUoKMjCV3Ivz7tGo0FJSQnWrl2LQ4cO4cyZM0hISIAkSejUqRO++OILjBo1Cg0aNLAoMa20PJLCcPLXX3/dNkbWBkdEpQDsH+YiPM/DZDLJdRzYDxuk5s2b3xbewf5mA6H8W2lVYa+zGnhvvvkmfv31Vxw5cgSOjo7QaDTIyMhAnz590LJlS/z000+wt7eH2WyWPbj9+/dHTEwM5s2bZ2E6JACF+UZcungBgiCB9PYo40WENwmBr4/PI40DojuYwUtKSpCXl4esrCwkJyfj0qVLSE9PR1ZWFtLT03Hz5k0LErDkFgaDQbZKubq6wtnZGU5OTvIJSBaz5efnhy+++AIeHh5yIrbU1FRs2LABAwYMQNu2bS3SuiprpVduP8dxyM3NxbBhw2A0GtG7d29ERkaif//+CAkJkT/L4t50Op3FCU7ldXieR6tWrfDBBx9gzJgxVkkSVxVqhCDVQVUTQ0mOqk4AKjsyPz8fERERGDduHN577z2UlpbCwcEBly5dwqBBg9C8eXOsWrUKbm5uclj84MGD0aNHD7z11luWtvXqFqyowT5gpeN+/fVXpKWl4dy5c0hJSUFubi5yc3PlU4T29vYICAhAgwYN4O7ujiZNmqBJkybw9/eHi4sLHBwc4OnpCRcXFxgMBuh0unuG6isnKVsEcnJy8PPPP6NRo0ZIS0sDz/MYM2aM7DepKqMj2xVYuYnKY3o/EgMbj5SUFAwYMACbNm1CmzZtrB6kqHyQGrNiVVYkK1ul2Od4nqeMjOuUe+MG5eRkk9FopKysLEpNTaWioiK6g8grW71WrVpFvr6+9NdffxERyabk9PR0ioqKorCwMNq+fbv8vYiICPr666+JSCHbSlSuaYpEoiASL4pkEs3Eiw/mma8OmNy+f/9+cnZ2ppYtW1Lfvn1pxowZ9OWXX9LGjRvp6NGjlJaWZlHm4GHHhynulT3pStN8u3btZD1yyJAhsm5XVZ8oFWxBECx0per2xc6dOyk0NFTWh2xBQSeqQTPv/aLcCiLS15+/Rz7ODtS73yhKu55Ln3/+Oc19ZR7dvHmDiASSJJ5IkizIwohHRDR48GAaM2YMEZV3MiNPSUkJzZ07l9zc3Cg6OppiYmIoLCyMLl68SESVQkgkkq8vElG5Olp7tvfS0lK6cePGXcNamCJcVlZGpaWlVFZWdlcFVhRFys/Pp+TkZDp9+jTFx8fTyZMnKTc3V+7DypOPXe/8+fMUGBhIDg4O5OjoSB4eHpSQkCB/r8rxVCyKDwJGkPnz51OPHj3k9tRmqM/dUOsiliQQNFqCYMrGyCeGYOfx61i2eiU4qRh9e/eEt5c7JKkMGk4LkK48ck6xQ7Ot988//8SQIUOwfPlyPPHEE+B53iJy9Ny5c1izZg1ycnIwdepURERE3Cb/2hIqn91QRghU1gHMZjPS09ORkpKCCxcuIDExEdevX0d6ejquX78OnudRWFiIsrIy2cixYsUKjB07tsrwDSZCFRUV4cknn8TevXsBlB9Q27ZtGxo3bvzIwj7YvSdNmgRJkrBy5UqLEg/WhhV0EAlm3gw7gz3OH9+DmJgnUeoUjK3bfkXn1qEo5c2w03Plde8kLZQJxBnYYM2ePRuHDx/G/v37YTAYLKJKlSWOAdxRt7E2SKF/3ek9o9GI48eP4/jx4zh69CjOnz+PzMxMCIIAZ2dn+Pj4oFmzZggKCkJAQAD8/Pzg5OQEOzs76HQ6ODo6olWrVnKGmMr3UoalX7x4EYsWLYJer0e/fv2QmJiIwYMHywWLarL/WFtEUUTv3r0xYsQIvPLKKzajoLNG1q6IJQlkFk3E82VEQgHNnTyCANC8f31OEhEV8zzxklAu+dxh12bbb0ZGBjVq1Eguz8VkaSYbm81mKisrk8WvhxEFrAEm+ixdupR0Oh0FBQXRE088Qf/4xz/ol19+oYSEBMrLy6tG30t3fU/pn2K/r169SvPnz5f9SErd8mHFIHaPjIwMCgwMpN27d1s8ty2g1gkiSBKZhWIi4mnvnh10Yt8emjzmKYLWibbGHSMiolJRrKh9V/U1WK0/ovI4IB8fH0pKSiq/PnMEKpRSonIn1KJFi4iodkPZHwasnZmZmZSQkHBbvFrlzzKPeuWf+3XMsn4tKyujlJQUyswsj0vLy8uj7777jn7//Xf5XpYhKg/3fAcOHKDAwEDZmWhL41PrBOElIqJS2rwplkZNmkTFxSKtWTSfAFBk32GUmpVDPBGZJYl4qnoTUZ4jISIaMWIE9enTR54IbPdQFgP19PSkAwcOyN+vC6hqIjKDxJ0sUg97P0mSqKioiLp160ZvvfUWEZF8v5UrV9LKlSstduqHvR8RUWxsLDVp0sSmQkwYal0g13EAkR4tW/bA239/FbxQgqgnRiAx8S8s+vpjuNkboAOg4zhoUbWbonI9jq+//hpXr17FrFmzANxK5GA2mzFz5kwsW7YMmzZtQlRUlEVqIFuH0pMuiqLsP9Lr9Rb+jpp6HuVxBOaVB27pChMnToSXlxcWLVqE/Px82fnI2vagSE5Olmsx0gOcg3mUsMqZdI7TIqxpS/nfbq7Ot3/mnte4VesvICAAGzduxKhRo9C/f39MmDABN27cwLJly8BxHPbu3Yv27dvDbDbbjvJXrf6qvUQSVFGTnJ3HAW6dC2cHsM6ePYtFixZh7NixCAoKks+TP+gBr/j4eDkH76MO76kurGbSIZJkyxIpzqNXdyFiyQDatm2LvXv3wtfXF/Pnz8cvv/yCyZMn49ChQ2jfvv1dTtGpUEKSJOh0Ojl2joHtXmazGa1bt8bo0aOxbt06xMfHy/VSqgO2UwiCgOzsbDRv3tzaj14lrJbVhOM0YAvOw4gJyoI2gYGBcHNzw9/+9jc5oyLzJ9hS8UtbhnIclMkU2OLCsjoGBwdj+vTpWLFiBbKzszFgwAA5P/D9BKQygly9ehX5+fno0KHDbfe3BdiWU6AGUFZWhqSkJAC4/3Q3KmQwXcLe3h43b94EcPukZc5HV1dXzJgxA2lpaVi7dq1FVa17gY1NWloaNBrNbZWJbQX1ZuawQQkODsaVK1cA2FbJ5LoIs9kMoOpJy9IIabVaTJ06FQ4ODli6dClMJtM967mTIqI7Pj4e/v7+8PPzszkFHahHBGEIDw/H1atX5cwoysFQcW+wvvL29kZeXt5dP8tKSUiShGHDhqFZs2ZYuHChnPqUhf9UhlKkTk1NRWBgoByqb2uoNwRhYlS7du1QWFiI8+fPA3iw9DmPM9gkdXFxua+qwsr8wD169EB0dDSWLVuGCxcuQK/X3/F7bLySk5PRqlWru97Dmqg3BGG7RfPmzeHj44MdO3YAgLqDVBNKUbWgoADFxcX3pcNxHAez2YxWrVph0qRJ2LJlC44ePSqLYkqRi5lys7KyEB8fL5c5sMVxqjcEAW51fI8ePbBz504Atqf02TKUUcOenp7IyclBQUHBfX2X5U6WJAn+/v6YMWMGTpw4ga1bt8oJqxlJ2O+ffvoJzs7O6N+/v3x/W0O9IgjD4MGDceHCBcTHxz9UrYzHEWwV9/HxgclkgtFotHj9XmC+EkdHR8yaNQt5eXn46aefLPQVnU4HQRCwYsUKjBs3Dk5OTjCbzdWyNla3tMODol4RhPk8oqKi0KRJE6xfvx6AbW7dtgo2ScPCwmBnZydbBKsDZbmFCRMmwN3dHQsXLkRZWZn8mT179iArKwsTJkyQv3M/kMwiRLE8S79IgFkUIImmCmczE6drbrzrHUEEQYBer8fUqVOxYcMGFBcXyzE+Ku4NRhB/f39wHIcLFy4AqN4iww6tsSiHgQMHomPHjli4cCGuX78OAFi+fDmioqIQEhIii8b3dC5KBI1Wi6VLFqBPr354YfY/kFWQj7g9v+H1V99FcVEpyhOiVrS1Boa8XhEEuDXAw4cPR2lpqbyLqGLW/YNV/G3atKlcBu9B9QNm4YqIiMCQIUOwbt06xMXF4fjx4xg/fnz1LsYB0AJ9ekXh1NH9OHzkJHw8vODm5oROnTrCxdWxPBi1BjNx1CuCsKA+URTh5eWFsWPHYsmSJfJ7tSW31nWwEJ0WLVrg3LlzchzWg/Qdi0jmeR6hoaGYMWMGFi5cCJ1Oh5iYGIsKX/dcxDgOJbyIsBYdMXvWRJw7sRsnTl+AMV9E164RFR+p2SldrwgCWFpiZs6cidTUVKxfv17OdK7i7lA68SIiIpCamirXhHwQgrDr6fV6CIIAJycnGAwGuLu74/jx4xb5se51fQ4AVzG2I8c9CS1M+N/Xi5FvtkPDoACQVPOWsHpHEOCWLhIcHIxnnnkGH330kUXpaBV3h7JOZHFxMc6cOQPg4YwdrLx3Wloa9u3bh3fffReFhYX44YcfANxvqTWCTmOGiUS069wLT/bsip+WLYeLjw+0eg1EUUJNW4rrJUGAW7rInDlzcPPmTSxcuBAA7lmGTMUtgjRt2hQNGzZEXFwcgKoJcr9iK/vMoUOH5AKrgwcPRsOGDREbG4uCgoJ7HsAicBCJA5nMgMEbA4fHoEV4GDq2Lz9LotWxDB8sTPzh+6LeEoR5dn18fPDaa6/h448/xrVr11SL1n2A6WsuLi6IjIzE/v375XzHVZ0RqU4dwW3btqFLly5yKEv//v0RFRWFBQsWyOEpdzuhqIVOTmPbK2Yg/t+c5+Ft0DyylE71liBs8CRJwqxZs9CoUSO89dZb8nsq7g7WRwMGDEBSUhKSk5OrVNTz8/NRXFx814lJFSmFSktLcfLkScTExMjv8TyPdu3a4bnnnsOWLVtw+PBh2SNfFfQaCV99+jE2/LIVRWVadOnCzpE8mjGttwTRaDQWDqvPP/8cv/32G7Zs2SJv5SruDDbhu3btCq1Wi8OHDwOwVKYlScLEiRPx5ZdfAriz+Mo+zzL9R0ZGAoB8vl6SJPj4+OD555/HyZMnsXXrVtnbbnFNIgAi9sXtxey/z8P5S1fRokVTiIIZ3COayvWWIMAts68gCOjatSvmzJmD2bNnIzc394GOidZFsAld3V2TESQ0NFTOsKh8nSE5ORnZ2dn3bAMAHDhwAA0bNkRoaCiAW95ztjPZ29vjhRdeQH5+PlasWCG/bkkSDZatXI6Vy5dgcPQTsLPTQZIeXQxXvSaIEqIo4u2330ZQUJCc/QSo/+JW5VID1fmeIJSHdPTv3x9HjhzBzZs3LWLbsrKykJeXJydcuNu1AGD//v2IiCj3V1TebViGRUEQMH78eHh7e+Obb76ByWSS78lxHCQRcHZxR+9ukXB3dYQg8iBOg0eVpv+xIAg7WajRaLBw4ULs378fn376qewbqW87iSRJsqI7b948fPDBB7LizQ4xVYcsAwcOlNOfArfqfezduxfFxcXo1q0bgKp9EDzPQ6PR4MaNG7hw4QJ69+59x/vodDpZSR80aBC6du2K2NhYZGVlyUGQnFYHTltRJEmrhb3WAL1e++jKWFAtJ46zJpSFPZ2cnORMgSwRWn2BsizB+++/TwBo9uzZcjZKVtbuXgnaWGK30tJSCg8PpylTplj0V3R0NHXv3p2I6J41Jvfu3UshISH3XaCT3SMlJYU++ugjudSFsrRCbWRgfKwIQnSrxsh7771Hfn5+crrLmkilaSuonHb1m2++IY1GQz179qRjx47JnyEiixSllaGsQfjee+/RtGnT5P7buXMncRxHK1eulPvvbv390UcfUa9eveR73quvWX1GovK6JV999ZVcNlpZGbc69S4fBI8dQZSr67hx46hDhw5yAuiysjKbygv7sFCW5j548CB17NiRANDgwYPp119/LS9WRLfIwvN8lUVwKhOuqKiIOnbsSH369LHI/1sV2HVGjhxJr7zyChHRbeW/7wblLrZgwQLasGGD/F5iYiLl5+erBKlJsMGWJIlKSkqoe/fuNGDAADKZTBYTqr5AWVyotLSUVqxYQW3btiWNRkMNGzaksWPH0rJly+TyzMp+kqvcViS1ZpP1ueeeIzc3N1nsUdZEr3wNovLk1+Hh4bRmzRr589UBy7VMRPTTTz/RunXr6OTJk9SgQQNasGCB/JyPAo8dQRg52GqYlZVFbdq0oZEjR9ZYUmZbAiO9coKbzWY6evQovfvuuxQREUF+fn7k5eVFERER9Pbbb9Phw4ctSq+xvsrPz6cpU6aQTqejjRs3yte6U3+x148cOWJR5au6k5klKmff27RpEwUEBBAA8vPzk6tgPQqSPHYEUYINfHp6OjVv3pxGjhxpUXOvvugkDFXVKJckia5fv07r16+nWbNmUZs2bUiv11OrVq3oiy++oBs3bhAR0V9//UXdu3cnvV4v7wT3quPBJuxXX30ll1d70Gz07DuZmZnUs2dPAkAODg4EgKKjo6msrOyOutTD4LEmiCiK8q5x+fJlCg8Pp7Fjx8qv2VIhl5oE2xVMJtNtz2gymWjbtm00evRosre3p1atWtFrr71GXl5eFBYWRocPH5b75n6tYBMmTKDnn3+eiB584WHfuXTpEs2ZM4d69OhBLi4uxJXHmND7778vt6sm8VgTRKmcEhElJydT8+bNacSIEXJ1WVZQsr7tJpUrRbHVV7nCnzt3jsaNG0eOjo700ksvyQV17meSs/dLSkqoQ4cO9MMPP9z3d+/nusXFxXT06FH69NNPqX///tSwYUPas2eP/Gw1hVqvUWirYHXxrly5gmHDhsHPzw9r1qyBp6cnTCYT9Hr9Y5HjlxRORI1GA5PJhMzMTAQFBcmO1ftJ6crqGSYkJGD48OH47bffEB4e/tB1Dln7lNcQBAGnTp0Cz/Po3r17jSYLrP8jfp9gAYzBwcHYvn07ioqKMHDgQKSlpcHOzq7eedvvBGV2dlZtNjg4GERUZQkJqjguWxmMZEeOHIGHhweaNGkiX78m2kdE4HkeZrMZWq0WXbp0Qbdu3Wq8QJJKkAqwilWiKCIgIABbt26Fr68v+vbtizNnzsjnqh+2mlJdAcdxcgUrNunYxGShLCwimq3mpDhfTooAxXbt2skh7DVZDUun08kR22xcaroGjEqQCrDIXxb96+7ujo0bN2LgwIGIjo7Gli1b5Dihx+VsOyMFS/rGiMB2DL1ej6ysLKxcuRJFRUUWB6d0Oh3Kysrw119/yfFXNb2wsDFTtrHGQY+xkl4VmINM6U3++OOPycPDgz7++GPZYVWffCX32yfKZ+Z5njZv3kzNmjWjBg0aWPg4mNHj2LFjFBoaSufPn5evU9egEqQKMEuOMhZp8+bN5O3tTU8++SRlZ2cT0S0nWU3UDL/fNjFfBpuIjMjMkcZCQth7ytfvVhJaSYKqnIvsO7m5ubRo0SLq0qULaTQaevLJJ+mvv/6yaB8zk//nP/+h7t2731aDvS5BFbGqgLLaK1B+duHJJ5/E0aNHUVhYiIiICMTFxckFQamS3P2owK7PREEWHs7+Zr/ZD6uGy15nR5DZST3lD8tLxfQEVilYo9GgrKwMf/zxB1544QW0bdsWzz//PDw9PbF9+3Zs2LAB4eHhcr+x9gHA7t270bt3bznU3haTU98Lqpn3LmATnw0w00/eeecdxMbG4vXXX8drr70mHw+9n/SZDwrWjpKSEqSlpSE9PR2ZmZlITU1FamoqSkpKAABeXl4IDw+Hj48PXF1dZfk8KCgIgYGBcHR0vOe9RFFEWloajh07hl27duHAgQO4evUqgoKC0L9/f0ycOFE+A8IIptfr5X7SaDRITU1F3759sWzZMvTu3VuuRlXXoBLkHmCrNjthx1bXrVu34oUXXkBYWBgWLlyIJk2ayJaU+y1kWR0w/8O2bdswYcIECIIAZ2dneHt7w87OTr4vq3FeUFAAnudRUlIi1xN0c3NDp06d0KFDB/j5+YHjOKSlpSE3Nxcmk0muSZicnIzU1FSYTCY0bdoUUVFRGDlyJCIjI28jGDsZyCY/I0JsbCwWL16MgwcPwmAw1NlCRipBqgm2q2i1Wly7dg1z5sxBXFwc/v3vf2P69OkW/oKatKowgqSnp+PYsWMIDg6Gv7+/TBBl+xhBCgsLUVhYCJPJhBs3biA+Ph6HDx/GpUuXUFJSIifSkyQJvr6+cHFxgYeHB4KCgtC+fXt07twZLVq0kK/PJnl+fj7WrFkDg8GAKVOmyCmBmBWL4zgMHDgQHTp0wEcffVRndw/20KqSXk1UVmC///578vf3p5iYGEpMTCQislB4awLMYPCwxgBJkig/P5+uXLlC169fp7y8PCooKLhnkF9xcTEdOnSIXnvtNWrSpAkBoBkzZsgKvdJQkZCQQMHBwfLhrLps8VN3kIcAU2p1Oh1SU1Mxd+5c7N+/H2+88QZefvll2Nvbyz4T5mR72JWUKulFbFf57rvv0LRpU3Ts2BFBQUGy/vGgEAQB165dQ0JCAnbu3Ik9e/bg0qVL8PLywuDBgzFp0iR07dpVTsTH2qXT6TBv3jz8+eef2L17d50VrRhUgjwkmFeZeZ1Xr16NN998E97e3vjPf/6Dfv36yWIXANnyVVP31mg0uHTpEqZNm4azZ88CABwdHdGsWTOEhIQgNDQUoaGh8PHxga+vLxwcHGSdhSr0q7KyMuTl5SEnJweZmZk4deoUTp06hcuXL8NkMiEwMBA9e/bE0KFD0bNnT7i7u9+xLdeuXUO3bt3w8ccfY8KECXVbvIJKkIeGsu4eEUGv1yMnJwcffPABVqxYgYEDB+Jf//oXmjVrJodEMEX/YVdXtmoz821GRgYSExNx5swZXLhwAYmJiTAajcjLy0NhYaFsEgZumWTZd0tLS6HT6eDm5gY/Pz+0b98eERER6Ny5M5o3bw57+1tTRBnWwa7DkoO/88472Lp1Kw4fPixbttQdRIUMZY7YI0eO4B//+Afi4+Mxbdo0zJkzB76+vhb+BiYGPegkIoWVrSqUlpbCaDQiKysLhYWFyMvLQ3Z2NsrKyqDRaODi4gIfHx94eXnBzc0NXl5e8PT0vE08YwuA0jrH7s2MFlevXkWvXr3w3nvvYfLkyXKwo0oQFRZgk4ntEps2bcL777+Pmzdv4v/9v/+HmTNnwtnZuUYtXkqRSUm+B5mcSnN1Zadp5Xuy+2m1WkyfPh1JSUkW2eBr2txd21AJ8ojAdgkW5GcymbBy5Up89tln4DgOL774IiZNmgRXV1eLkPGH3VEAy6DAytnYlas/M8sqyUAKb/09awZWEFyv12Pr1q2YMmUKNm/ejG7duslZGev6GRqVII8YlSdgQUEBvv/+e8TGxkKn02HmzJlyqk2lMv/IolNrEMw3c+3aNfTq1QtTpkzBu+++W+cVcyVUgtQSlOZZjUaDwsJCLFu2DLGxsSgpKcGYMWMwY8YMNGnSxCLhtK2KKGx3LCoqQkxMDDw9PbFp0yYAlmXw6jpUgtQiSBHUyESYwsJC/PLLL1i8eDEuXbqEmJgYTJ8+HREREfIkY+EczJdS25OPFNEDyvYXFRVh5MiRKCwsxJYtW+Dh4VGvdg9AJYjVUFn04nkef/zxB7777jscOXIErVq1wrPPPouhQ4fCw8MDwK3AQBaZW5ttVUYr63Q6XL16VfZzrF+/HgEBAQ993twWoRLEyqjKN3Lu3DksX74cmzZtgk6nw9ChQzF+/Hi0a9futuOtteFnUO5gHMdhw4YNmD17NiIiIrB48eJ6uXMwqASxAShXaEYUAMjLy8PGjRuxatUqXLx4EWFhYRg7diyGDBmCwMBA+fuCIMiiV1U1A+9FIKXvRmm6VepAHMfhxIkT+PDDD3H48GG8+uqrmDt3ruxotFVd6WGhEsQGUdnkS0RISEjAjz/+iN9//x0FBQXo2LEjhg8fjn79+sHf31/+Ljv8pJzYlSeukjxKEzMAC4IC5eTbv38/YmNjERcXh+7du+Odd95Bhw4dHlnhTFuCShAbBpvsSpNvSUkJDh48iA0bNiAuLg6iKCIiIgJPPPEE+vTpg6CgIItrKE/yKcPR74bCwkIkJCQgLi4OW7duRWpqKrp27YqXXnoJffr0kdtWX3cNJVSC2DCUzrzK1i+gfCIfOHAAv/76K/bt24fS0lK0bNkS/fv3x4ABA9C0aVOLsyJKmM1mFBUVobCwENnZ2bh48SIuX76MhIQEJCYmoqioCCEhIejfvz9GjRqF5s2bA7AUvYCHz3Nl61AJUgdRWQQDysly4sQJbN26Ffv370dubi68vLzg6upqYYFizkhBEFBaWori4mIAgKurKzw9PdG6dWt0794dHTp0QOPGjWVxS+mXeZygEqQO405h9IIg4Pz58zh27Bhyc3PB87wcOGhvbw+tVgtvb2/4+/vD09MTHh4e8Pf3v+04LVVkL9RoNPeVbrQ+QiVIHQcTwyor2g8i+igVdqXeovz9uEElSD2E0lR7N6gkuDcez32znqNyfioVD47HS+NSoaKaUAmiQsVdoBJEhYq7QAdABEAAHo8KMSpUVAM6AHoAHABVo1OhohJ0AA6g3Mxb/8smqVChQoWKmsP/B7Zv6kshDu2PAAAAJXRFWHRkYXRlOmNyZWF0ZQAyMDE5LTA0LTA0VDAwOjQ3OjIxLTA3OjAwFJLj/AAAACV0RVh0ZGF0ZTptb2RpZnkAMjAxOS0wNC0wNFQwMDo0NzoyMS0wNzowMGXPW0AAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from IPython.display import Image\n", "Image(filename='Data/coordinatesystem_neuromag.png')\n", "#downloaded from url='http://www.fieldtriptoolbox.org/_media/faq/coordinatesystem_neuromag.png?w=200&tok=e89b5d')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "MNI units are milimeters. The 3d image data are transformed by the associated affine transformation such that\n", " to match an averaged brain template, called MNI152 (is it s an average brain taken from MRI's of 152 healthy individuals). \n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The code below is the first attempt to use Plotly for selection and brain slice plotting:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import plotly.graph_objs as go\n", "import plotly.tools as tls" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import copy\n", "from matplotlib import cm as mpl_cm\n", "import matplotlib as mpl\n", "from nilearn import (plotting, _utils)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def mpl_to_plotly(cmap, pl_entries):\n", " # convert a mpl colormap to a Plotly colorscale with pl_entries (colors)\n", " h=1.0/(pl_entries-1)\n", " pl_colorscale=[]\n", " for k in range(pl_entries):\n", " C=list(map(np.uint8, np.array(cmap(k*h)[:3])*255))\n", " pl_colorscale.append([round(k*h,2), f'rgb({C[0]}, {C[1]}, {C[2]})'])\n", " return pl_colorscale" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def colorscale(cmap, values, threshold=None, symmetric_cmap=True,\n", " vmax=None, vmin=None):\n", " \"\"\"This function modifies nilearn.plotting.js_plotting_utils.colorscale:\n", " https://github.com/nilearn/nilearn/blob/master/nilearn/plotting/js_plotting_utils.py#L192\n", " It defines a Plotly colorscale or a mixed Plotly colorscale from a given nilearn or matplotlib colormap,\n", " and extracts the color range \"\"\"\n", " \n", " cmap = mpl_cm.get_cmap(cmap)\n", " abs_values = np.abs(values)\n", " \n", " if not symmetric_cmap and (values.min() < 0):\n", " warnings.warn('you have specified symmetric_cmap=False '\n", " 'but the map contains negative values; '\n", " 'setting symmetric_cmap to True')\n", " symmetric_cmap = True\n", " if symmetric_cmap and vmin is not None:\n", " warnings.warn('vmin cannot be chosen when cmap is symmetric')\n", " vmin = None\n", " if threshold is not None:\n", " if vmin is not None:\n", " warnings.warn('choosing both vmin and a threshold is not allowed; '\n", " 'setting vmin to 0')\n", " vmin = 0\n", " if vmax is None:\n", " vmax = abs_values.max()\n", " if symmetric_cmap:\n", " vmin = - vmax\n", " if vmin is None:\n", " vmin = values.min()\n", " \n", " #Define a mixed colorscale consisting from 3 parts, with gray at center, and and two diverging cmaps outside\n", " L=11\n", " abs_threshold = None\n", " if threshold is not None:\n", " abs_threshold = _utils.param_validation.check_threshold(threshold, values, _utils.extmath.fast_abs_percentile)\n", " norm = mpl.colors.Normalize(vmin=vmin, vmax=vmax)\n", " ca = norm(-abs_threshold)\n", " cb = norm(abs_threshold)\n", " h1 = ca/L\n", " dl = [k*h1 for k in range(L)]\n", " h2 = (cb-ca) / (L-1)\n", " dc = [ca+k*h2 for k in range(L)]\n", " h3 = (1-cb-h2) / (L-1)\n", " dr = [cb+h2+k*h3 for k in range(L)]\n", " d = dl+dc+dr\n", " cmaplist = [cmap(t)[:3] for t in d]\n", " for k in range(L): \n", " cmaplist[k+L] = mpl_cm.gray(k*0.1)[:3]\n", " \n", " pl_colorscale = []\n", " \n", " for k, t in enumerate(d):\n", " c = list(map(np.uint8, np.array(cmaplist[k])*255))\n", " pl_colorscale.append([round(t,3), f'rgb({c[0]}, {c[1]}, {c[2]})']) \n", " else:\n", " pl_colorscale = mpl_to_plotly(cmap, L)\n", " return {\n", " 'colorscale': pl_colorscale, 'vmin': vmin, 'vmax': vmax, \n", " 'abs_threshold': abs_threshold}" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def pl_view_img(stat_map_img, bg_img='MNI152', \n", " threshold=1e-6,\n", " cmap=plotting.cm.cold_hot,\n", " symmetric_cmap=True,\n", " dim='auto',\n", " vmax=None,\n", " vmin=None,\n", " resampling_interpolation='continuous',\n", " **kwargs):\n", " \"\"\"\n", " Reads and processes the background image (MNI152) and a stats image, to be able to\n", " perform an interactive Plotly view of the coresponding statistical map;\n", " THIS IS a MODIFIED version of the nilearn view_img function:\n", " https://github.com/nilearn/nilearn/blob/master/nilearn/plotting/html_stat_map.py#L332\n", " \n", " Returns\n", " -------\n", " color_info, background, statistical and mask NIfTI images,\n", " the affine transformation associated to the processed stat image,\n", " \n", " \n", " \"\"\" \n", " \n", " mask_img, stat_map_img, data, threshold = plotting.html_stat_map._mask_stat_map(stat_map_img, threshold)\n", " color_info = colorscale(cmap, data.ravel(), threshold=threshold, \n", " symmetric_cmap=symmetric_cmap, vmax=vmax,\n", " vmin=vmin)\n", " \n", " bg_img, bg_min, bg_max, black_bg = plotting.html_stat_map._load_bg_img(stat_map_img, bg_img, dim)\n", " \n", " stat_map_img, mask_img = plotting.html_stat_map._resample_stat_map(stat_map_img, bg_img, mask_img,\n", " resampling_interpolation)\n", " \n", " affine_transf = stat_map_img.affine\n", " \n", " return [color_info, \n", " _utils.niimg._safe_get_data(bg_img, ensure_finite=True), \n", " stat_map_img, \n", " _utils.niimg._safe_get_data(mask_img, ensure_finite=True), \n", " affine_transf]\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def get_plotly_heatmaps(bg_img, stat_map_img, mask_img, affine_transf, color_info,\n", " cut_coords= None, threshold=3):\n", " \n", " # Extract data from the three NIfTI images needed to define the Plotly heatmap traces\n", " # representing the brain slices\n", " \n", " # bg_img: a numpy array storing the background image data\n", " # stat_map_img: a NIfTI image, NOT its data yet\n", " # mask_img: numpy array; mask image data\n", " # affine_transf: numpy array of shape (4,4) as the affine transformation of the stat_map_img\n", " # cut_coords: 3-list giving the MNI coordinates where the brain cut is performed\n", " # color_info: dict for plotly colorscale, and vmin, vmax values to be colormapped\n", " \n", " # Returns \n", " # the Plotly heatmap traces \n", " # and the min-max values of MNI coordinates in each direction\n", " \n", " cut_slices = plotting.html_stat_map._get_cut_slices(stat_map_img, cut_coords, threshold)\n", " \n", " stat_map_img = _utils.niimg._safe_get_data(stat_map_img, ensure_finite=True) \n", " \n", " \n", " #min-max MNI-coordinates corresponding to min-max voxel indices i, j, k\n", " xMNI_min, yMNI_min, zMNI_min = affine_transf[:, 3][:3]\n", " imax, jmax, kmax = stat_map_img.shape \n", " xMNI_max, yMNI_max, zMNI_max, one = np.dot(affine_transf, [imax-1, jmax-1, kmax-1, 1])\n", " \n", " pl_colorscale=color_info['colorscale']\n", " vmin=color_info['vmin']\n", " vmax=color_info['vmax']\n", " abs_threshold=color_info['abs_threshold']\n", " \n", " islice, jslice, kslice = np.array(cut_slices-1, int)# voxel indices corresponding to cut_slices\n", " \n", " # Mix the backgraound and statistical image values according to mask image:\n", " a, b = -abs_threshold, abs_threshold\n", " vmin_bg, vmax_bg = bg_img.min(), bg_img.max()\n", " final_img = copy.deepcopy(stat_map_img) #define the image that will be plotted\n", " alpha = (b-a) / (vmax_bg-vmin_bg) \n", " new_bg = a + alpha * (bg_img-vmin_bg) # map bg_img vals to [a,b]\n", " I, J, K = np.where(mask_img==1)\n", " final_img[I, J, K] = new_bg[I, J, K]\n", " \n", " # Define 2d arrays storing the image values in the three slices:\n", " xsts = final_img[islice, :, :]\n", " ysts = final_img[:, jslice, :]\n", " zsts = final_img[:, :, kslice]\n", " \n", " # Trace definition: The second char in the x, y, z names\n", " # represents the slice name\n", " #x-slice:\n", " xx = np.linspace(yMNI_min, yMNI_max, jmax)\n", " yx = np.linspace(zMNI_min, zMNI_max, kmax)\n", " textx = [[f'y: {xx[j]}
z: {yx[i]}
val: '+'{:.2f}'.format(xsts.T[i,j]) for j in range(jmax)] \n", " for i in range(kmax)]\n", "\n", " slicex = dict(type='heatmap',\n", " zsmooth='best',\n", " x=xx,\n", " y=yx,\n", " z=xsts.T,\n", " text=textx, \n", " colorscale=pl_colorscale, \n", " colorbar=dict(thickness=20, ticklen=4, tick0=-7, dtick=2),\n", " showscale=False, \n", " hoverinfo='text', \n", " zmin=vmin,\n", " zmax=vmax)\n", " # y-slice\n", " xy = np.linspace(xMNI_min, xMNI_max, imax)\n", " yy = np.linspace(zMNI_min, zMNI_max, kmax)\n", " texty = [[f'x: {xy[j]}
z: {yy[i]}
val: '+'{:.2f}'.format(ysts.T[i,j]) for j in range(imax)]\n", " for i in range(kmax)]\n", "\n", " slicey = dict(type='heatmap',\n", " x=xy,\n", " y=yy, \n", " z=ysts.T,\n", " text=texty, \n", " zsmooth='best', \n", " colorscale=pl_colorscale, \n", " colorbar=dict(thickness=20, ticklen=4, tick0=-7, dtick=2),\n", " showscale=False, \n", " hoverinfo='text', \n", " zmin=vmin,\n", " zmax=vmax)\n", " # z-slice\n", " xz = np.linspace(xMNI_min, xMNI_max, imax)\n", " yz = np.linspace(yMNI_min, yMNI_max, jmax)\n", " textz = [[f'x: {xz[j]}
y: {yz[i]}
val: '+'{:.2f}'.format(zsts.T[i,j]) for j in range(imax)]\n", " for i in range(jmax)]\n", "\n", " slicez = dict(type='heatmap',\n", " zsmooth='best',\n", " x=xz,\n", " y=yz, \n", " z=zsts.T,\n", " text=textz, \n", " colorscale=pl_colorscale, \n", " colorbar=dict(thickness=20, ticklen=4, tick0=np.ceil(vmin), dtick=2),\n", " hoverinfo='text',\n", " zmin=vmin,\n", " zmax=vmax)\n", " return [slicex, slicey, slicez, xMNI_min, yMNI_min, zMNI_min, xMNI_max, yMNI_max, zMNI_max] " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Interactive view of brain slices through a point given in MNI-coordinates" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We define a subplot of three heatmaps representing the orthogonal brain slices through a given point in MNI-coordinates. \n", "A mouse click on a slice (subplot cell) retrieves the MNI coordinates that define a new cut point in the brain: \n", " a coordinate is the coordinate of the clicked slice (let us say y in the case of \n", "the y-slice), and the coordinates of the clicked point (x and z in our example), that complement the first one.\n", "\n", "A callback function executed `on_click` retrieves a new cut point and generates the new slices through that point. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "stat_img = 'Data/image_10426.nii.gz'# from nilearn datasets\n", "\n", "color_info, bg_img, stat_map_img, mask_img, affine_transf = pl_view_img(stat_img, threshold=3)\n", "\n", "cut_coords = [13, -27, 50]# Important!!! define cut_coords as a list (not as a tuple!!) because it is\n", " # updated during the interactivity\n", " \n", "slicex, slicey, slicez, xMNI_min, yMNI_min, zMNI_min, xMNI_max, yMNI_max, zMNI_max = \\\n", " get_plotly_heatmaps(bg_img, stat_map_img, mask_img, affine_transf, color_info,\n", " cut_coords= cut_coords, threshold=3)\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The affine transformation, associated to the processed `stat_map_img` is listed below:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "affine_transf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It consists in a scaling of factor 2, and a translation that places the voxel of index (0,0,0) to the point of MNI coordinates\n", "$(x=-90, y=-126, z=-72)$.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define the Plotly subplots:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fig= tls.make_subplots(subplot_titles=(f'Slice x={cut_coords[0]}', f'Slice y={cut_coords[1]}',\n", " f'Slice z={cut_coords[2]}'), \n", " rows=1,\n", " cols=3,\n", " horizontal_spacing=0.065,\n", " )" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fw = go.FigureWidget(fig) # convert fig to a FigureWidget " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fw.add_trace(slicex, 1, 1)\n", "fw.add_trace(slicey, 1, 2)\n", "fw.add_trace(slicez, 1, 3)#HOW TO STOP THIS automaticaly json display?" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with fw.batch_update():\n", " fw.layout.update(title='Orthogonal brain slices through a point',\n", " width=1000, height=425)\n", " fw.layout.xaxis1.update(range=[yMNI_min, yMNI_max], title='y')\n", " fw.layout.yaxis1.update(range=[zMNI_min, zMNI_max], title='z')\n", " fw.layout.xaxis2.update(range=[xMNI_min, xMNI_max], title='x')\n", " fw.layout.yaxis2.update(range=[zMNI_min, zMNI_max], title='z')\n", " fw.layout.xaxis3.update(range=[xMNI_min, xMNI_max], title='x')\n", " fw.layout.yaxis3.update(range=[yMNI_min, yMNI_max], title='y')\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now let us define a callback function to be executed when a point in one of the three slices is clicked." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "slice_name=['x', 'y', 'z']\n", "\n", "def clicked_trace(trace, points, state):\n", " global cut_coords, slice_name\n", " if points.trace_index in [0, 1, 2]:\n", " idx = points.trace_index\n", " if len(points.xs) != 0:\n", " xclicked, yclicked = points.xs[0], points.ys[0]\n", " if idx == 0:\n", " cut_coords[1:] = [xclicked, yclicked]\n", " elif idx == 1:\n", " cut_coords[0], cut_coords[2] = xclicked, yclicked \n", " else: # idx == 2:\n", " cut_coords[:2] = [xclicked, yclicked ] \n", " tr = get_plotly_heatmaps(bg_img, stat_map_img, mask_img, affine_transf, color_info,\n", " cut_coords= cut_coords)\n", " with fw.batch_update():\n", " fw.update(data=[tr[0], tr[1], tr[2]])\n", " for k in range(3):\n", " fw.layout.annotations[k].text=f'Slice {slice_name[k]}={cut_coords[k]}'\n", "\n", "for trace in fw.data: \n", " trace.on_click(clicked_trace)\n", "fw " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let display the `nilearn` version of the above slices:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "\n", "view = plotting.view_img(stat_img, cut_coords=cut_coords, threshold=3)\n", "view" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The Plotly plot displays the range of axes, and provides the possibility to select a particular point \n", "in a slice, by hovering an area of interest, and choosing the point of desired coordinates,\n", "while in the original nilearn plot a blind point selection is performed, whose coordinates are seen a posteriori." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To realize how Plotly cut point selection works, before running this notebook, see the gif file:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%%html\n", "" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.1" } }, "nbformat": 4, "nbformat_minor": 2 }