{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"I love doing data analyses with pandas, numpy, sci-py etc., but I often need to run repeated measures ANOVAs, which are not implemented in any major python libraries. Python Psychologist shows how to do repeated measures ANOVAs yourself in python, but I find using a widley distributed implementation comforting... \n",
"\n",
"In this post I show how to execute a repeated measures ANOVAs using the rpy2 library, which allows us to move data between python and R, and execute R commands from python. I use rpy2 to load the car library and run the ANOVA. \n",
"\n",
"I will show how to run a one-way repeated measures ANOVA and a two-way repeated measures ANOVA. "
]
},
{
"cell_type": "code",
"execution_count": 86,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"#first import the libraries I always use. \n",
"import numpy as np, scipy.stats, pandas as pd\n",
"\n",
"import matplotlib as mpl\n",
"import matplotlib.pyplot as plt\n",
"import pylab as pl\n",
"%matplotlib inline\n",
"pd.options.display.mpl_style = 'default'\n",
"plt.style.use('ggplot')\n",
"mpl.rcParams['font.family'] = ['Bitstream Vera Sans']\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Below I use the random library to generate some fake data. I seed the random number generator with a one so that this analysis can be replicated. \n",
"\n",
"I will generated 3 conditions which represent 3 levels of a single variable.\n",
"\n",
"The data are generated from a gaussian distribution. The second condition has a higher mean than the other two conditions."
]
},
{
"cell_type": "code",
"execution_count": 92,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXUAAAECCAYAAAAB2kexAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFKVJREFUeJzt3V9MW/fdx/GP/zxJGv4c97gwFmgeaEFquxYWzURDbIGm\nVfWonQbzRSfBpoVxsa3TpFrTlioSbSS0RVUecFBHlotO6y53U3w19c6gCW5MxNSsSqqHdd2EImLh\nf4G0pLHNcxHFa1qCj4MN5df36w6fP/46P/Xdw+GEuDY2NjYEADCCe7cHAACUD1EHAIMQdQAwCFEH\nAIMQdQAwCFEHAIN4i+1w9epVnTt3rvD1tWvX9P3vf1/Hjh1TOBzWysqK6urqFAqFVFVVJUmamppS\nNBqV2+3W0NCQOjo6KvcJAAAFrlKeU8/n8/rpT3+q3/72t3rnnXdUU1Ojvr4+RSIR3bhxQ4ODg1pa\nWtLExITOnDmjZDKp0dFRTUxMyO3mmwIAqLSSSnvp0iU1NDTooYce0vz8vHp6eiRJvb29isVikqRY\nLKbu7m55vV7V19eroaFBi4uL5Z8cAPA5JUV9dnZW3d3dkqRMJiOfzydJsixLmUxGkpRKpeT3+wvH\n+P1+JZPJcs0LANiC46hns1ldvHhRXV1dn9vmcrm2PLbYdgBAeTiO+sLCgh555BHV1tZKun11nk6n\nJd2+OrcsS5Jk27YSiUThuEQiIdu2yzkzAOAeij79csenb71IUiAQ0PT0tPr7+zUzM6POzs7C6xMT\nE/rOd76jZDKp5eVltba23vO8V69e3cb4X2w1NTVaXV3d7TFwn1i/vcv0tTt06NA9tzmK+vr6ui5d\nuqSf/OQnhdf6+/sVDocVjUYLjzRKUlNTk7q6uhQKheTxeDQ8PMztFwDYISU90lgJXKnji4r127tM\nX7utrtR5eBwADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADOL4\ntzQC25G8uaH42q3dHqMknnRauVxut8dwrL76v2Tv55fnfdkRdeyI+NotnXzng90ew2iv/88jsvfv\n2+0xsMu4/QIABiHqAGAQog4ABiHqAGAQog4ABiHqAGAQog4ABiHqAGAQog4ABiHqAGAQR78m4MaN\nG7pw4YKWlpYkSS+99JK++tWvKhwOa2VlRXV1dQqFQqqqqpIkTU1NKRqNyu12a2hoSB0dHZX7BACA\nAkdR/+Mf/6gjR47ol7/8pXK5nG7evKm3335b7e3t6uvrUyQSUSQS0eDgoJaWljQ3N6fx8XElk0mN\njo5qYmJCbjffFABApRUt7UcffaQrV67o+PHjkiSPx6ODBw9qfn5ePT09kqTe3l7FYjFJUiwWU3d3\nt7xer+rr69XQ0KDFxcUKfgQAwB1Fr9Tj8bhqa2t1/vx5/etf/1JLS4tOnDihTCYjn88nSbIsS5lM\nRpKUSqXU1tZWON7v9yuZTG57UH51687g17cCe1vRqOdyOf3zn//Uj3/8Y7W2tuqtt95SJBK5ax+X\na+sIbLW9pqbG0aD/l07zq1t3wP++0Kb/fshX9vN60umynxN383g8jv97Mt2+ffu+tH8WRaPu9/tl\n27ZaW1slSd/85jc1NTUln8+ndDotn8+nVColy7IkSbZtK5FIFI5PJBKybfue519dXXU06F674t2r\ncrmc4zUp9byorEqt3V5UU1Nj9J/FVv/DKnpP3efz6aGHHtLVq1clSe+++64efvhhfeMb39D09LQk\naWZmRp2dnZKkQCCg2dlZZbNZxeNxLS8vF/6HAACoLEdPvwwNDemNN95QNpvVV77yFb300kvK5/MK\nh8OKRqOFRxolqampSV1dXQqFQvJ4PBoeHi56ewbAFxc/z9oZ5fp5lqOoNzc368yZM597fWRkZNP9\ng8GggsHg9iYD8IXAP0W4M8r1zxHy8DgAGISoA4BBiDoAGISoA4BBiDoAGISoA4BBiDoAGISoA4BB\niDoAGISoA4BBiDoAGISoA4BBiDoAGISoA4BBiDoAGISoA4BBiDoAGISoA4BBiDoAGISoA4BBiDoA\nGISoA4BBiDoAGMTrZKef//zneuCBB+R2u+XxeHTmzBmtra0pHA5rZWVFdXV1CoVCqqqqkiRNTU0p\nGo3K7XZraGhIHR0dFf0QAIDbHEVdkk6fPq3q6urC15FIRO3t7err61MkElEkEtHg4KCWlpY0Nzen\n8fFxJZNJjY6OamJiQm433xQAQKU5Lu3GxsZdX8/Pz6unp0eS1Nvbq1gsJkmKxWLq7u6W1+tVfX29\nGhoatLi4WMaRAQD34uhK3eVyaXR0VG63W88++6yeffZZZTIZ+Xw+SZJlWcpkMpKkVCqltra2wrF+\nv1/JZLICowMAPstR1EdHR/Xggw/q+vXrGh0dVWNj413bXS7Xlsdvtb2mpsbJCPKk0472w/Z4PB7H\na1LSeVm/imPt9rZyrZ+jqD/44IOSpNraWh09elSLi4uyLEvpdFo+n0+pVEqWZUmSbNtWIpEoHJtI\nJGTb9j3Pvbq66mjQXC7naD9sTy6Xc7wmpZ4XlcXa7W2lrN9W8S96T/3mzZv6+OOPJUnr6+t69913\ndfjwYQUCAU1PT0uSZmZm1NnZKUkKBAKanZ1VNptVPB7X8vKyWltbHQ0KANieolfqmUxGZ8+elSTl\n83l961vfUkdHhx599FGFw2FFo9HCI42S1NTUpK6uLoVCIXk8Hg0PDxe9PQMAKI+iUa+vry9E/dOq\nq6s1MjKy6THBYFDBYHD70wEASsLD4wBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh\n6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBg\nEKIOAAYh6gBgEK+TnfL5vF555RXZtq1XXnlFa2trCofDWllZUV1dnUKhkKqqqiRJU1NTikajcrvd\nGhoaUkdHR0U/AADgPxxdqf/lL39RU1OTXC6XJCkSiai9vV0TExN68sknFYlEJElLS0uam5vT+Pi4\nTp06pTfffFP5fL5y0wMA7lI06olEQgsLCzp+/Lg2NjYkSfPz8+rp6ZEk9fb2KhaLSZJisZi6u7vl\n9XpVX1+vhoYGLS4uVnB8AMCnFY36n/70J/3gBz+Q2/2fXTOZjHw+nyTJsixlMhlJUiqVkt/vL+zn\n9/uVTCbLPTMA4B62vKd+8eJF1dbWqqWlRe+9996m+9y5JXMvxbbX1NQUGfE2TzrtaD9sj8fjcbwm\nJZ2X9as41m5vK9f6bRn1999/XxcvXtTCwoJu3bqljz/+WG+88YYsy1I6nZbP51MqlZJlWZIk27aV\nSCQKxycSCdm2veUAq6urjgbN5XKO9sP25HI5x2tS6nlRWazd3lbK+m0V/y1vvwwMDOj3v/+9Jicn\n9fLLL+trX/uafvGLXygQCGh6elqSNDMzo87OTklSIBDQ7Oysstms4vG4lpeX1dra6vAjAQC2y9Ej\njXfcuZXS39+vcDisaDRaeKRRkpqamtTV1aVQKCSPx6Ph4eGit18AAOXjOOpPPPGEnnjiCUlSdXW1\nRkZGNt0vGAwqGAyWZzoAQEn4G6UAYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoA\nYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoAYBCiDgAGIeoAYBCi\nDgAGIeoAYBDvVhs/+eQTnT59Wrdu3VI2m1VnZ6cGBga0tramcDislZUV1dXVKRQKqaqqSpI0NTWl\naDQqt9utoaEhdXR07MgHAQAUifq+ffv02muvaf/+/crlcnr11Vd15coVzc/Pq729XX19fYpEIopE\nIhocHNTS0pLm5uY0Pj6uZDKp0dFRTUxMyO3mGwIA2AlFa7t//35JUjabVT6fV1VVlebn59XT0yNJ\n6u3tVSwWkyTFYjF1d3fL6/Wqvr5eDQ0NWlxcrOD4AIBP2/JKXZLy+bxOnjypa9eu6bnnntPDDz+s\nTCYjn88nSbIsS5lMRpKUSqXU1tZWONbv9yuZTFZodADAZxWNutvt1tmzZ/XRRx/pN7/5jf7+97/f\ntd3lcm15fLHtNTU1DsaUPOm0o/2wPR6Px/GalHRe1q/iWLu9rVzrVzTqdxw8eFBHjhzRBx98IMuy\nlE6n5fP5lEqlZFmWJMm2bSUSicIxiURCtm1ved7V1VVH75/L5ZyOim3I5XKO16TU86KyWLu9rZT1\n2yr+W95Tv379um7cuCHp9pMwly5dUktLiwKBgKanpyVJMzMz6uzslCQFAgHNzs4qm80qHo9reXlZ\nra2tjoYEAGzfllfq6XRak5OTyufz2tjY0LFjx/TUU0+ppaVF4XBY0Wi08EijJDU1Namrq0uhUEge\nj0fDw8NFb78AAMpny6gfPnxYr7/++uder66u1sjIyKbHBINBBYPB8kwHACgJD5ADgEGIOgAYhKgD\ngEGIOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEGI\nOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEGIOgAYhKgDgEG8xXZYWVnR5OSkMpmMXC6XnnnmGT3/\n/PNaW1tTOBzWysqK6urqFAqFVFVVJUmamppSNBqV2+3W0NCQOjo6Kv5BAAAOou71evWjH/1Izc3N\nWl9f18mTJ9Xe3q7p6Wm1t7err69PkUhEkUhEg4ODWlpa0tzcnMbHx5VMJjU6OqqJiQm53XxTAACV\nVrS0Pp9Pzc3NkqQDBw6osbFRyWRS8/Pz6unpkST19vYqFotJkmKxmLq7u+X1elVfX6+GhgYtLi5W\n7hMAAApKunyOx+P68MMP1dbWpkwmI5/PJ0myLEuZTEaSlEql5Pf7C8f4/X4lk8kyjgwAuJeit1/u\nWF9f19jYmE6cOKEHHnjgrm0ul2vLY7faXlNT4+j9Pem0o/2wPR6Px/GalHRe1q/iWLu9rVzr5yjq\n2WxWY2NjOnbsmI4ePSrp9tV5Op2Wz+dTKpWSZVmSJNu2lUgkCscmEgnZtn3Pc6+urjoaNJfLOdoP\n25PL5RyvSannRWWxdntbKeu3VfyL3n7Z2NjQhQsX1NjYqBdeeKHweiAQ0PT0tCRpZmZGnZ2dhddn\nZ2eVzWYVj8e1vLys1tZWR4MCALan6JX6+++/r7/+9a86fPiwfv3rX0uSBgYG1N/fr3A4rGg0Wnik\nUZKamprU1dWlUCgkj8ej4eHhordnAADlUTTqjz32mP785z9vum1kZGTT14PBoILB4PYmAwCUjIfH\nAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAg\nRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAgRB0ADELUAcAg3mI7nD9/XgsLC6qt\nrdXY2JgkaW1tTeFwWCsrK6qrq1MoFFJVVZUkaWpqStFoVG63W0NDQ+ro6KjsJwAAFBS9Un/66ad1\n6tSpu16LRCJqb2/XxMSEnnzySUUiEUnS0tKS5ubmND4+rlOnTunNN99UPp+vzOQAgM8pGvXHH3+8\ncBV+x/z8vHp6eiRJvb29isVikqRYLKbu7m55vV7V19eroaFBi4uLFRgbALCZ+7qnnslk5PP5JEmW\nZSmTyUiSUqmU/H5/YT+/369kMlmGMQEAThS9p16My+Xa1vaamhpH7+NJpx3PhPvn8Xgcr0lJ52X9\nKo6129vKtX73FXXLspROp+Xz+ZRKpWRZliTJtm0lEonCfolEQrZtb3mu1dVVR++Zy+XuZ1SUKJfL\nOV6TUs+LymLt9rZS1m+r+N/X7ZdAIKDp6WlJ0szMjDo7Owuvz87OKpvNKh6Pa3l5Wa2trffzFgCA\n+1D0Sv3cuXO6fPmyrl+/rp/97Gd68cUX1d/fr3A4rGg0WnikUZKamprU1dWlUCgkj8ej4eHhordf\nAADlUzTqL7/88qavj4yMbPp6MBhUMBjc3lQAgPvC3ygFAIMQdQAwCFEHAIMQdQAwCFEHAIMQdQAw\nCFEHAIMQdQAwCFEHAIMQdQAwCFEHAIMQdQAwCFEHAIMQdQAwCFEHAIMQdQAwCFEHAIMQdQAwCFEH\nAIMQdQAwCFEHAIMQdQAwiLcSJ/3b3/6mt956S/l8XsePH1d/f38l3gYA8Bllv1LP5/P6wx/+oFOn\nTml8fFyzs7NaWloq99sAADZR9qgvLi6qoaFB9fX18nq96u7u1vz8fLnfBgCwibJHPZlMyu/3F762\nbVvJZLLcbwMA2AQ/KAUAg5T9B6W2bSuRSBS+TiQSsm37nvsfOnTI0XkPHZJiTzVvdzzsEtZv72Lt\n9payX6k/+uijWl5eVjweVzab1dzcnAKBQLnfBgCwCdfGxsZGuU+6sLBw1yON3/ve98r9FgCATVQk\n6gCA3cEPSgHAIEQdAAxSkV8T8GV3/vx5LSwsqLa2VmNjY7s9DkqwsrKiyclJZTIZuVwuPfPMM3r+\n+ed3eyw49Mknn+j06dO6deuWstmsOjs7NTAwsNtj7SjuqVfA5cuXdeDAAf3ud78j6ntMOp1WOp1W\nc3Oz1tfXdfLkSf3qV79SU1PTbo8Gh27evKn9+/crl8vp1Vdf1Q9/+EM99thjuz3WjuH2SwU8/vjj\nqqqq2u0xcB98Pp+am5slSQcOHFBjY6NSqdTuDoWS7N+/X5KUzWaVz+dVXV29yxPtLG6/APcQj8f1\n4Ycfqq2tbbdHQQny+bxOnjypa9eu6bnnnvvSfZfFlTqwifX1dY2Pj+vEiRM6cODAbo+DErjdbp09\ne1YXLlzQ5cuX9d577+32SDuKqAOfkc1mNTY2pm9/+9s6evTobo+D+3Tw4EEdOXJE//jHP3Z7lB1F\n1IFP2djY0IULF9TY2KgXXnhht8dBia5fv64bN25Iuv0kzKVLl9TS0rLLU+0snn6pgHPnzuny5cta\nXV2VZVl68cUX9fTTT+/2WHDgypUreu2113T48GG5XC5J0sDAgL7+9a/v8mRw4t///rcmJyeVz+e1\nsbGhY8eO6bvf/e5uj7WjiDoAGITbLwBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAYh6gBgEKIOAAb5\nf3RZLDXn/fnfAAAAAElFTkSuQmCC\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"import random\n",
"\n",
"random.seed(1) #seed random number generator\n",
"cond_1 = [random.gauss(600,30) for x in range(30)] #condition 1 has a mean of 600 and standard deviation of 30\n",
"cond_2 = [random.gauss(650,30) for x in range(30)] #u=650 and sd=30\n",
"cond_3 = [random.gauss(600,30) for x in range(30)] #u=600 and sd=30\n",
"\n",
"plt.bar(np.arange(1,4),[np.mean(cond_1),np.mean(cond_2),np.mean(cond_3)],align='center') #plot data\n",
"plt.xticks([1,2,3]);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, I load rpy2 for ipython. I am doing these analyses with ipython in a jupyter notebook (highly recommended). \n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"%load_ext rpy2.ipython"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here's how to run the ANOVA. Note that this is a one-way anova with 3 levels of the factor. "
]
},
{
"cell_type": "code",
"execution_count": 97,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"Type III Repeated Measures MANOVA Tests:\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: (Intercept) \n",
"\n",
" Response transformation matrix:\n",
" (Intercept)\n",
"cond_1 1\n",
"cond_2 1\n",
"cond_3 1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" (Intercept)\n",
"(Intercept) 102473990\n",
"\n",
"Sum of squares and products for error:\n",
" (Intercept)\n",
"(Intercept) 78712.7\n",
"\n",
"Multivariate Tests: (Intercept)\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.9992 37754.33 1 29 < 2.22e-16 ***\n",
"Wilks 1 0.0008 37754.33 1 29 < 2.22e-16 ***\n",
"Hotelling-Lawley 1 1301.8736 37754.33 1 29 < 2.22e-16 ***\n",
"Roy 1 1301.8736 37754.33 1 29 < 2.22e-16 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: Factor \n",
"\n",
" Response transformation matrix:\n",
" Factor1 Factor2\n",
"cond_1 1 0\n",
"cond_2 0 1\n",
"cond_3 -1 -1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" Factor1 Factor2\n",
"Factor1 3679.584 19750.87\n",
"Factor2 19750.870 106016.58\n",
"\n",
"Sum of squares and products for error:\n",
" Factor1 Factor2\n",
"Factor1 40463.19 27139.59\n",
"Factor2 27139.59 51733.12\n",
"\n",
"Multivariate Tests: Factor\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.7152596 35.16759 2 28 2.303e-08 ***\n",
"Wilks 1 0.2847404 35.16759 2 28 2.303e-08 ***\n",
"Hotelling-Lawley 1 2.5119704 35.16759 2 28 2.303e-08 ***\n",
"Roy 1 2.5119704 35.16759 2 28 2.303e-08 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"Univariate Type III Repeated-Measures ANOVA Assuming Sphericity\n",
"\n",
" SS num Df Error SS den Df F Pr(>F) \n",
"(Intercept) 34157997 1 26238 29 37754.334 < 2.2e-16 ***\n",
"Factor 59964 2 43371 58 40.094 1.163e-11 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"\n",
"Mauchly Tests for Sphericity\n",
"\n",
" Test statistic p-value\n",
"Factor 0.96168 0.57866\n",
"\n",
"\n",
"Greenhouse-Geisser and Huynh-Feldt Corrections\n",
" for Departure from Sphericity\n",
"\n",
" GG eps Pr(>F[GG]) \n",
"Factor 0.96309 2.595e-11 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
" HF eps Pr(>F[HF])\n",
"Factor 1.03025 1.163294e-11\n",
"\n"
]
}
],
"source": [
"#pop the data into R\n",
"%Rpush cond_1 cond_2 cond_3 \n",
"\n",
"#label the conditions\n",
"%R Factor <- c('Cond1','Cond2','Cond3')\n",
"#create a vector of conditions\n",
"%R idata <- data.frame(Factor) \n",
"\n",
"#combine data into single matrix\n",
"%R Bind <- cbind(cond_1,cond_2,cond_3) \n",
"#generate linear model\n",
"%R model <- lm(Bind~1)\n",
"\n",
"#load the car library. note this library must be installed.\n",
"%R library(car) \n",
"#run anova\n",
"%R analysis <- Anova(model,idata=idata,idesign=~Factor,type=\"III\") \n",
"#create anova summary table\n",
"%R anova_sum = summary(analysis) \n",
"\n",
"#move the data from R to python\n",
"%Rpull anova_sum \n",
"print anova_sum"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The ANOVA table isn't pretty, but it works. As you can see, the ANOVA was wildly significant. \n",
"\n",
"Next, I generate data for a two-way (2x3) repeated measures ANOVA. Condition A is the same data as above. Condition B has a different pattern (2 is lower than 1 and 3), which should produce an interaction. "
]
},
{
"cell_type": "code",
"execution_count": 98,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXUAAAECCAYAAAAB2kexAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGEZJREFUeJzt3V1sXFe99/HvvDRO4peZztSuiZ0oSW2dNoBNVDvFchu7\naVWdtkljXJ6CEhAxFgKKkDpCkCjCTSUL2qrEE6tNyAWFVoILlAiPQEKlvRhbyBZiHJmnoSTlMaWA\nlcfx8bzFTu3E83Iucjqnpc7MONneblZ+nzvvPXuv/97L88ue5ZU1jmw2m0VERIzgXOkCRETEOgp1\nERGDKNRFRAyiUBcRMYhCXUTEIAp1ERGDuAu94Ny5cxw5ciT38/nz5/nCF77A9u3bCQaDTE9PU1lZ\nSSAQoLS0FICBgQHC4TBOp5Ouri4aGxuX7wpERCTHsZR56plMhm984xv88Ic/5LXXXqO8vJzdu3cT\nCoW4ePEie/fuZWJigv7+fp599llisRi9vb309/fjdOpDgYjIcltS0p4+fZrq6mpuu+02RkdHaWtr\nA6C9vZ1IJAJAJBKhtbUVt9tNVVUV1dXVjI+PW1+5iIh8xJJCfXh4mNbWVgCSySRerxcAj8dDMpkE\nIB6P4/f7c8f4/X5isZhV9YqISB5Fh3oqleLUqVO0tLR8ZJ/D4ch7bKH9IiJijaJDfWxsjM2bN1NR\nUQFceTpPJBLAladzj8cDgM/nIxqN5o6LRqP4fD4raxYRkasoOPvlfR8cegFoampicHCQjo4OhoaG\naG5uzm3v7+9n586dxGIxJicnqauru+p5z507dx3lf7yVl5czMzOz0mXINVL/3bhM77t169ZddV9R\noT4/P8/p06f5+te/ntvW0dFBMBgkHA7npjQC1NbW0tLSQiAQwOVy0d3dreEXERGbLGlK43LQk7p8\nXKn/blym912+J3VNHhcRMYhCXUTEIAp1ERGDKNRFRAyiUBcRMYhCXUTEIAp1ERGDKNRFRAxS9DIB\nIiI3ioXJc7gmJ2xrz3Hb7aQ8H481rhTqImKczH9Ncvm5/ba1t+rA8/AxCXUNv4iIGERP6nLDiF3K\nMjW7YFt7Nek5yl22NSdiCYW63DCmZhfY/9o7trX3o0frKfcq1eXGouEXERGDKNRFRAyiUBcRMYhC\nXUTEIAp1ERGDKNRFRAyiUBcRMYhCXUTEIAp1ERGDKNRFRAxS1DIBFy9e5Pjx40xMXFnK8sknn+QT\nn/gEwWCQ6elpKisrCQQClJaWAjAwMEA4HMbpdNLV1UVjY+PyXYGIiOQUFeo/+9nP2Lp1K9/5zndI\np9NcunSJX/3qVzQ0NLB7925CoRChUIi9e/cyMTHByMgIfX19xGIxent76e/vx+nUhwIRkeVWMGnf\ne+89zp49y44dOwBwuVysXbuW0dFR2traAGhvbycSiQAQiURobW3F7XZTVVVFdXU14+Pjy3gJIiLy\nvoJP6lNTU1RUVHDs2DH+8Y9/sGnTJvbt20cymcTr9QLg8XhIJpMAxONx6uvrc8f7/X5isdgylb80\nWrpVRExXMNTT6TR///vf+epXv0pdXR2vvPIKoVDoQ69xOBx5z5Fvf3l5eZGlXr//l0jYunTr4V3/\nwbqaCtvaM50rkbC1PYfTaevvp1hnwZk/k6zmcrlZ+zH5XSkY6n6/H5/PR11dHQCf/exnGRgYwOv1\nkkgk8Hq9xONxPB4PAD6fj2g0mjs+Go3i8139a55mZmau9xqKlk6nbWsLIJvJ2Hp9plP/SbFWZbK2\ntpdOp2z9Xcn3sFEw1L1eL7fddhvnzp1j3bp1vPnmm6xfv57169czODhIR0cHQ0NDNDc3A9DU1ER/\nfz87d+4kFosxOTmZ+wfhZlOzkMD1tynb2vs4ffmtiKyMoma/dHV18eKLL5JKpbj99tt58sknyWQy\nBINBwuFwbkojQG1tLS0tLQQCAVwuF93d3QWHZ0zljk1x+Uc355ffinyQ3X/PujNl76e6j5OiQn3j\nxo08++yzH9ne09Oz6Os7Ozvp7Oy8vspExBh2fxXhiaYsN+scBU0eFxExiEJdRMQgCnUREYMo1EVE\nDKJQFxExiEJdRMQgCnUREYMo1EVEDKJQFxExiEJdRMQgRS0TIHIzsnNBNi3GJlZRqItchZ0Lsmkx\nNrGKhl9ERAyiUBcRMYhCXUTEIAp1ERGDKNRFRAyiUBcRMYhCXUTEIAp1ERGDKNRFRAyiUBcRMYhC\nXUTEIEWt/fKtb32LNWvW4HQ6cblcPPvss8zOzhIMBpmenqayspJAIEBpaSkAAwMDhMNhnE4nXV1d\nNDY2LutFiIjIFUUv6PXMM89QVlaW+zkUCtHQ0MDu3bsJhUKEQiH27t3LxMQEIyMj9PX1EYvF6O3t\npb+/H6dTHwpERJZb0UmbzWY/9PPo6ChtbW0AtLe3E4lEAIhEIrS2tuJ2u6mqqqK6uprx8XELSxYR\nkasp6knd4XDQ29uL0+nkwQcf5MEHHySZTOL1egHweDwkk0kA4vE49fX1uWP9fj+xWGwZShcRkX9X\nVKj39vZy6623cuHCBXp7e6mpqfnQfofDkff4fPvLy8uLKcESrkTCtrYAyH9bLOdyuVlr4/20m8n9\np76z2E383isq1G+99VYAKioq2LZtG+Pj43g8HhKJBF6vl3g8jsfjAcDn8xGNRnPHRqNRfL6rL/4/\nMzNzPfUvSTqdtq0tALKFX2KldDpl6/20m8n9p76zmOHvvXwPwwXH1C9dusTc3BwA8/PzvPnmm2zY\nsIGmpiYGBwcBGBoaorm5GYCmpiaGh4dJpVJMTU0xOTlJXV2dBZchIiKFFHxSTyaTvPDCCwBkMhnu\nvfdeGhsbueOOOwgGg4TD4dyURoDa2lpaWloIBAK4XC66u7sLDs+IiIg1CoZ6VVVVLtQ/qKysjJ6e\nnkWP6ezspLOz8/qrExGRJdHkcRERgyjURUQMolAXETGIQl1ExCAKdRERgyjURUQMolAXETGIQl1E\nxCAKdRERgyjURUQMolAXETGIQl1ExCAKdRERgyjURUQMolAXETGIQl1ExCAKdRERgyjURUQMolAX\nETGIQl1ExCAKdRERgyjURUQMolAXETGIu5gXZTIZDhw4gM/n48CBA8zOzhIMBpmenqayspJAIEBp\naSkAAwMDhMNhnE4nXV1dNDY2LusFiIjI/yrqSf23v/0ttbW1OBwOAEKhEA0NDfT39/OpT32KUCgE\nwMTEBCMjI/T19XHw4EF+8pOfkMlklq96ERH5kIKhHo1GGRsbY8eOHWSzWQBGR0dpa2sDoL29nUgk\nAkAkEqG1tRW3201VVRXV1dWMj48vY/kiIvJBBUP91Vdf5Utf+hJO5/++NJlM4vV6AfB4PCSTSQDi\n8Th+vz/3Or/fTywWs7pmERG5irxj6qdOnaKiooJNmzbx1ltvLfqa94dkrqbQ/vLy8gIlWseVSNjW\nFgD5L91yLpebtTbeT7uZ3H/qO4vdxO+9vKH+9ttvc+rUKcbGxlhYWGBubo4XX3wRj8dDIpHA6/US\nj8fxeDwA+Hw+otFo7vhoNIrP58tbwMzMjAWXUZx0Om1bWwBk7W0unU7Zej/tZnL/qe8sZvh7L9/D\ncN7hlz179vDjH/+Yo0eP8tRTT/HJT36Sb3/72zQ1NTE4OAjA0NAQzc3NADQ1NTE8PEwqlWJqaorJ\nyUnq6uqsuxIREcmrqCmN73t/KKWjo4NgMEg4HM5NaQSora2lpaWFQCCAy+Wiu7u74PCLiIhYp+hQ\n37JlC1u2bAGgrKyMnp6eRV/X2dlJZ2enNdWJiMiS6H+UiogYRKEuImIQhbqIiEEU6iIiBlGoi4gY\nRKEuImIQhbqIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEEU6iIi\nBlGoi4gYRKEuImIQhbqIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEHc+XZevnyZZ555hoWFBVKpFM3N\nzezZs4fZ2VmCwSDT09NUVlYSCAQoLS0FYGBggHA4jNPppKuri8bGRlsuRERECoT6qlWrOHToECUl\nJaTTaZ5++mnOnj3L6OgoDQ0N7N69m1AoRCgUYu/evUxMTDAyMkJfXx+xWIze3l76+/txOvWBQETE\nDgXTtqSkBIBUKkUmk6G0tJTR0VHa2toAaG9vJxKJABCJRGhtbcXtdlNVVUV1dTXj4+PLWL6IiHxQ\n3id1gEwmw/79+zl//jwPPfQQ69evJ5lM4vV6AfB4PCSTSQDi8Tj19fW5Y/1+P7FYbJlKFxGRf1cw\n1J1OJy+88ALvvfceP/jBD/jzn//8of0OhyPv8YX2l5eXF1GmNVyJhG1tAZD/0i3ncrlZa+P9tJvJ\n/ae+s9hN/N4rGOrvW7t2LVu3buWdd97B4/GQSCTwer3E43E8Hg8APp+PaDSaOyYajeLz+fKed2Zm\n5hpLX7p0Om1bWwBk7W0unU7Zej/tZnL/qe8sZvh7L9/DcN4x9QsXLnDx4kXgykyY06dPs2nTJpqa\nmhgcHARgaGiI5uZmAJqamhgeHiaVSjE1NcXk5CR1dXUWXYaIiBSS90k9kUhw9OhRMpkM2WyW7du3\n8+lPf5pNmzYRDAYJh8O5KY0AtbW1tLS0EAgEcLlcdHd3Fxx+ERER6+QN9Q0bNvD8889/ZHtZWRk9\nPT2LHtPZ2UlnZ6c11YmIyJJoArmIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEEU6iIiBlGoi4gYRKEu\nImIQhbqIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEEU6iIiBlGoi4gYRKEuImIQhbqIiEEU6iIiBlGo\ni4gYRKEuImIQhbqIiEEU6iIiBlGoi4gYxF3oBdPT0xw9epRkMonD4eCBBx7gkUceYXZ2lmAwyPT0\nNJWVlQQCAUpLSwEYGBggHA7jdDrp6uqisbFx2S9ERESKCHW3281XvvIVNm7cyPz8PPv376ehoYHB\nwUEaGhrYvXs3oVCIUCjE3r17mZiYYGRkhL6+PmKxGL29vfT39+N06kOBiMhyK5i0Xq+XjRs3ArB6\n9WpqamqIxWKMjo7S1tYGQHt7O5FIBIBIJEJraytut5uqqiqqq6sZHx9fvisQEZGcJT0+T01N8e67\n71JfX08ymcTr9QLg8XhIJpMAxONx/H5/7hi/308sFrOwZBERuZqCwy/vm5+f5/Dhw+zbt481a9Z8\naJ/D4ch7bL795eXlxZZw3VyJhG1tAZD/tljO5XKz1sb7aTeT+099Z7Gb+L1XVKinUikOHz7M9u3b\n2bZtG3Dl6TyRSOD1eonH43g8HgB8Ph/RaDR3bDQaxefzXfXcMzMz11P/kqTTadvaAiBrb3PpdMrW\n+2k3k/tPfWcxw997+R6GCw6/ZLNZjh8/Tk1NDY8++mhue1NTE4ODgwAMDQ3R3Nyc2z48PEwqlWJq\naorJyUnq6uqu8xJERKQYBZ/U3377bX7/+9+zYcMGvve97wGwZ88eOjo6CAaDhMPh3JRGgNraWlpa\nWggEArhcLrq7uwsOz4iIiDUKhvqdd97JL3/5y0X39fT0LLq9s7OTzs7O66tMRESWTJPHRUQMolAX\nETGIQl1ExCAKdRERgyjURUQMolAXETGIQl1ExCAKdRERgyjURUQMolAXETGIQl1ExCBFr6cuInKj\nKKtdz6pD/eC0ZzFBx4bNsKbU0nOmUinm5uaWfJxCXUSM0rzBi+tWH7OlZfY1msqAxeupr1mzhltu\nuYWFhYUlHafhFxExyv9p/ARzly+vdBnXbW5ujtWrVy/5OIW6iBilxHVzf3+DQl1EzHKTfymPQl1E\nxCAKdRERg2j2i4gYL3Ypy9Ts0maRLEVV2S34Soof9vn85z/PmTNnGBsbY9WqVZbWolAXEeNNzS6w\n/7V3lu38z//nZnwlxYXzv/71L8bGxqipqeH1119n586dltai4RcRERudPHmS++67j8cff5wTJ05Y\nfn6FuoiIjU6ePMljjz3Grl27GBoaYnp62tLzK9RFRGzyxz/+kcnJSR566CE2b95MfX09AwMDlrZR\ncEz92LFjjI2NUVFRweHDhwGYnZ0lGAwyPT1NZWUlgUCA0tIr6x4MDAwQDodxOp10dXXR2NhoacEi\nIjeqEydOsH37dsrKrixhsHPnTk6cOMHXvvY1y9ooGOr3338/Dz/8MC+99FJuWygUoqGhgd27dxMK\nhQiFQuzdu5eJiQlGRkbo6+sjFovR29tLf38/Tqc+EIjIzW1ubo7f/OY3ZDIZtm7dCsDly5dJJpP8\n5S9/YcuWLZa0UzBt77rrrtxT+PtGR0dpa2sDoL29nUgkAkAkEqG1tRW3201VVRXV1dWMj49bUqiI\nyI3sd7/7HS6Xi8HBQd544w3eeOMNBgcHueeeezh58qRl7VzTlMZkMonX6wXA4/GQTCYBiMfj1NfX\n517n9/uJxWIWlCkicu2qym7h+f/cvKznL+TkyZN88YtfZN26dR/avm/fPg4dOsT3v/99S0Y1rnue\nuqPAOguF9peXl19vCUVzJRK2tQWAzUtQuFxu1tp4P+1mcv+p75aXr8RR9Dzy5fLzn/980e27du1i\n165di+5zuVxLzshrCnWPx0MikcDr9RKPx/F4PAD4fD6i0WjuddFoFJ/Pl/dcMxavQZxPOp22rS0A\nsvY2l06nbL2fdjO5/9R3sph0Or3o70W+oL+mZ/2mpiYGBwcBGBoaorm5Obd9eHiYVCrF1NQUk5OT\n1NXVXUsTIiJyDQo+qR85coQzZ85w4cIFvvnNb/LEE0/Q0dFBMBgkHA7npjQC1NbW0tLSQiAQwOVy\n0d3dXXD4RURErFMw1J966qlFt/f09Cy6vbOzk87OzuurSkRErokmkIuIGEShLiJiEIW6iIhBFOoi\nIgZRqIuIGETffCQixnMnY2Snzy/b+R233U7Kk/8/Wt5zzz1MT0/jcrm45ZZbuPvuu3nuuec+smzA\n9VKoi4jxstPnufzc/mU7/6oDz0OBUHc4HLz66qvce++9XLp0iYMHD9LT08PLL79saS0afhERsVlJ\nSQmPPPIIf/3rXy0/t0JdRMQm2eyVBYXm5ub49a9/zd133215Gxp+ERGxQTabpbu7G7fbzXvvvYff\n7+cXv/iF5e0o1EVEbOBwOPjpT3/KvffeSzab5bXXXuPxxx9ncHCQyspKy9rR8IuIiM0cDgcPP/ww\nLpcr981xVtGTuoiITd4fU89ms7z++uskk8kPfVucFRTqImI8x223X5l2uIznL8a+fftwuVw4HA7W\nr19Pf3+/Ql1EZKlSHl/BeeTL7Q9/+IMt7WhMXUTEIAp1ERGDKNRFRAyiUBcRs/zPDJOblUJdRIxy\nKa1QFxExxon/+/9Zs2rVSpdx3dasWcP8/PySj9OURhExSuSfCdLxMsri/wVOhy1tOjZshjWllp4z\nlUqxsLCw5OOWJdT/9Kc/8corr5DJZNixYwcdHR3L0YyIyKJmJ/6F60fLt376v1t14HnSd9xlW3v5\nWD78kslkePnllzl48CB9fX0MDw8zMTFhdTMiIrIIy0N9fHyc6upqqqqqcLvdtLa2Mjo6anUzIiKy\nCMtDPRaL4ff7cz/7fD5isZjVzYiIyCI0+0VExCCW/6HU5/MRjUZzP0ejUXy+qy+kY/U3aeezbh1E\nPr3RtvYAuF9DT1ZR/9241Hf2sfxJ/Y477mBycpKpqSlSqRQjIyM0NTVZ3YyIiCzCkc1a/39qx8bG\nPjSl8XOf+5zVTYiIyCKWJdRFRGRl6A+lIiIGUaiLiBhEa78sg2PHjjE2NkZFRQWHDx9e6XJkCaan\npzl69CjJZBKHw8EDDzzAI488stJlSZEuX77MM888w8LCAqlUiubmZvbs2bPSZdlKY+rL4MyZM6xe\nvZqXXnpJoX6DSSQSJBIJNm7cyPz8PPv37+e73/0utbW1K12aFOnSpUuUlJSQTqd5+umn+fKXv8yd\nd9650mXZRsMvy+Cuu+6itNTaFdvEHl6vl40bNwKwevVqampqiMfjK1uULElJSQlwZZXDTCZDWVnZ\nCldkLw2/iFzF1NQU7777LvX19StdiixBJpNh//79nD9/noceeuim+5SlJ3WRRczPz9PX18e+fftY\nvXr1SpcjS+B0OnnhhRc4fvw4Z86c4a233lrpkmylUBf5N6lUisOHD3Pfffexbdu2lS5HrtHatWvZ\nunUrf/vb31a6FFsp1EU+IJvNcvz4cWpqanj00UdXuhxZogsXLnDx4kXgykyY06dPs2nTphWuyl6a\n/bIMjhw5wpkzZ5iZmcHj8fDEE09w//33r3RZUoSzZ89y6NAhNmzYgMNx5avQ9uzZw2c+85kVrkyK\n8c9//pOjR4+SyWTIZrNs376dxx57bKXLspVCXUTEIBp+ERExiEJdRMQgCnUREYMo1EVEDKJQFxEx\niEJdRMQgCnUREYMo1EVEDPLfIJhI1n3rDJ8AAAAASUVORK5CYII=\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"random.seed(1)\n",
"\n",
"cond_1a = [random.gauss(600,30) for x in range(30)] #u=600,sd=30\n",
"cond_2a = [random.gauss(650,30) for x in range(30)] #u=650,sd=30\n",
"cond_3a = [random.gauss(600,30) for x in range(30)] #u=600,sd=30\n",
"\n",
"cond_1b = [random.gauss(600,30) for x in range(30)] #u=600,sd=30\n",
"cond_2b = [random.gauss(550,30) for x in range(30)] #u=550,sd=30\n",
"cond_3b = [random.gauss(650,30) for x in range(30)] #u=650,sd=30\n",
"\n",
"width = 0.25\n",
"plt.bar(np.arange(1,4)-width,[np.mean(cond_1a),np.mean(cond_2a),np.mean(cond_3a)],width)\n",
"plt.bar(np.arange(1,4),[np.mean(cond_1b),np.mean(cond_2b),np.mean(cond_3b)],width,color=plt.rcParams['axes.color_cycle'][0])\n",
"plt.legend(['A','B'],loc=4)\n",
"plt.xticks([1,2,3]);"
]
},
{
"cell_type": "code",
"execution_count": 84,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"Type III Repeated Measures MANOVA Tests:\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: (Intercept) \n",
"\n",
" Response transformation matrix:\n",
" (Intercept)\n",
"cond_1a 1\n",
"cond_2a 1\n",
"cond_3a 1\n",
"cond_1b 1\n",
"cond_2b 1\n",
"cond_3b 1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" (Intercept)\n",
"(Intercept) 401981075\n",
"\n",
"Sum of squares and products for error:\n",
" (Intercept)\n",
"(Intercept) 185650.5\n",
"\n",
"Multivariate Tests: (Intercept)\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.9995 62792.47 1 29 < 2.22e-16 ***\n",
"Wilks 1 0.0005 62792.47 1 29 < 2.22e-16 ***\n",
"Hotelling-Lawley 1 2165.2575 62792.47 1 29 < 2.22e-16 ***\n",
"Roy 1 2165.2575 62792.47 1 29 < 2.22e-16 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: Factor1 \n",
"\n",
" Response transformation matrix:\n",
" Factor11\n",
"cond_1a 1\n",
"cond_2a 1\n",
"cond_3a 1\n",
"cond_1b -1\n",
"cond_2b -1\n",
"cond_3b -1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" Factor11\n",
"Factor11 38581.51\n",
"\n",
"Sum of squares and products for error:\n",
" Factor11\n",
"Factor11 142762.3\n",
"\n",
"Multivariate Tests: Factor1\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.2127533 7.837247 1 29 0.0090091 **\n",
"Wilks 1 0.7872467 7.837247 1 29 0.0090091 **\n",
"Hotelling-Lawley 1 0.2702499 7.837247 1 29 0.0090091 **\n",
"Roy 1 0.2702499 7.837247 1 29 0.0090091 **\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: Factor2 \n",
"\n",
" Response transformation matrix:\n",
" Factor21 Factor22\n",
"cond_1a 1 0\n",
"cond_2a 0 1\n",
"cond_3a -1 -1\n",
"cond_1b 1 0\n",
"cond_2b 0 1\n",
"cond_3b -1 -1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" Factor21 Factor22\n",
"Factor21 91480.01 77568.78\n",
"Factor22 77568.78 65773.02\n",
"\n",
"Sum of squares and products for error:\n",
" Factor21 Factor22\n",
"Factor21 90374.60 56539.06\n",
"Factor22 56539.06 87589.85\n",
"\n",
"Multivariate Tests: Factor2\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.5235423 15.38351 2 28 3.107e-05 ***\n",
"Wilks 1 0.4764577 15.38351 2 28 3.107e-05 ***\n",
"Hotelling-Lawley 1 1.0988223 15.38351 2 28 3.107e-05 ***\n",
"Roy 1 1.0988223 15.38351 2 28 3.107e-05 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"------------------------------------------\n",
" \n",
"Term: Factor1:Factor2 \n",
"\n",
" Response transformation matrix:\n",
" Factor11:Factor21 Factor11:Factor22\n",
"cond_1a 1 0\n",
"cond_2a 0 1\n",
"cond_3a -1 -1\n",
"cond_1b -1 0\n",
"cond_2b 0 -1\n",
"cond_3b 1 1\n",
"\n",
"Sum of squares and products for the hypothesis:\n",
" Factor11:Factor21 Factor11:Factor22\n",
"Factor11:Factor21 179585.9 384647\n",
"Factor11:Factor22 384647.0 823858\n",
"\n",
"Sum of squares and products for error:\n",
" Factor11:Factor21 Factor11:Factor22\n",
"Factor11:Factor21 92445.33 45639.49\n",
"Factor11:Factor22 45639.49 89940.37\n",
"\n",
"Multivariate Tests: Factor1:Factor2\n",
" Df test stat approx F num Df den Df Pr(>F) \n",
"Pillai 1 0.901764 128.5145 2 28 7.7941e-15 ***\n",
"Wilks 1 0.098236 128.5145 2 28 7.7941e-15 ***\n",
"Hotelling-Lawley 1 9.179605 128.5145 2 28 7.7941e-15 ***\n",
"Roy 1 9.179605 128.5145 2 28 7.7941e-15 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"Univariate Type III Repeated-Measures ANOVA Assuming Sphericity\n",
"\n",
" SS num Df Error SS den Df F Pr(>F) \n",
"(Intercept) 66996846 1 30942 29 62792.4662 < 2.2e-16 ***\n",
"Factor1 6430 1 23794 29 7.8372 0.009009 ** \n",
"Factor2 26561 2 40475 58 19.0310 4.42e-07 ***\n",
"Factor1:Factor2 206266 2 45582 58 131.2293 < 2.2e-16 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"\n",
"Mauchly Tests for Sphericity\n",
"\n",
" Test statistic p-value\n",
"Factor2 0.96023 0.56654\n",
"Factor1:Factor2 0.99975 0.99648\n",
"\n",
"\n",
"Greenhouse-Geisser and Huynh-Feldt Corrections\n",
" for Departure from Sphericity\n",
"\n",
" GG eps Pr(>F[GG]) \n",
"Factor2 0.96175 6.876e-07 ***\n",
"Factor1:Factor2 0.99975 < 2.2e-16 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
" HF eps Pr(>F[HF])\n",
"Factor2 1.028657 4.420005e-07\n",
"Factor1:Factor2 1.073774 2.965002e-22\n",
"\n"
]
}
],
"source": [
"%Rpush cond_1a cond_1b cond_2a cond_2b cond_3a cond_3b\n",
"\n",
"%R Factor1 <- c('A','A','A','B','B','B')\n",
"%R Factor2 <- c('Cond1','Cond2','Cond3','Cond1','Cond2','Cond3')\n",
"%R idata <- data.frame(Factor1, Factor2)\n",
"\n",
"#make sure the vectors appear in the same order as they appear in the dataframe\n",
"%R Bind <- cbind(cond_1a, cond_2a, cond_3a, cond_1b, cond_2b, cond_3b)\n",
"%R model <- lm(Bind~1)\n",
"\n",
"%R library(car)\n",
"%R analysis <- Anova(model, idata=idata, idesign=~Factor1*Factor2, type=\"III\")\n",
"%R anova_sum = summary(analysis)\n",
"%Rpull anova_sum\n",
"\n",
"print anova_sum"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Again, the anova table isn't too pretty. \n",
"\n",
"This obviously isn't the most exciting post in the world, but its a nice bit of code to have in your back pocket if you're doing experimental analyses in python."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"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.11"
}
},
"nbformat": 4,
"nbformat_minor": 0
}