face-detection.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import sqlite3\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "from matplotlib import pyplot as plt\n",
    "import seaborn as sns\n",
    "sns.set_style('white')\n",
    "%matplotlib inline\n",
    "\n",
    "from bgfunc import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dbfile = 'project_june10_705pm.db'\n",
    "conn = sqlite3.connect(dbfile)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "platform = 'ig' # ig = instagram, tw = twitter\n",
    "condition = 'ptsd' # depression, pregnancy, ptsd, cancer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "specs = analysis_specifications(platform, condition)\n",
    "\n",
    "platform_long = specs['plong'][platform]\n",
    "gb_types = specs['gb_types'][platform]\n",
    "fields = specs['fields'][platform] \n",
    "test_name = specs['test_name'][condition]\n",
    "test_cutoff = specs['test_cutoff'][condition]\n",
    "photos_rated = specs['photos_rated'][condition]\n",
    "has_test = specs['has_test'][condition]\n",
    "\n",
    "clfs = ['rf'] # lr = logistic regression, rf = random forests, svc = support vector\n",
    "periods = ['before']\n",
    "turn_points = ['from_diag']\n",
    "\n",
    "#################\n",
    "#\n",
    "# INSTAGRAM ONLY \n",
    "# \n",
    "#################\n",
    "additional_data = False # (instagram only) do we also get photo ratings? (This is proxy for limit_date_range in places)\n",
    "\n",
    "include_filter = True # (instagram only) include binary indicator of filter use? \n",
    "\n",
    "include_face_data = True # (instagram only) include has_face and face_ct? \n",
    "\n",
    "get_ratings_correlation = False # (instagram only) gets inter-rater correlation averages for each rated variable \n",
    "\n",
    "compare_ig_filters = False # (instagram only) compares use of filters between target and control pops\n",
    "\n",
    "#################\n",
    "#\n",
    "# TWITTER ONLY \n",
    "# \n",
    "#################\n",
    "populate_wordfeats_db = False # (twitter only) generates word features from reagan code \n",
    "\n",
    "#################\n",
    "#\n",
    "# ALL PLATFORMS \n",
    "# \n",
    "#################\n",
    "impose_test_cutoff = True # do we want to limit target pop based on testing cutoff (eg. cesd > 21)?\n",
    "\n",
    "limit_date_range = False # do we want to restrict posts by date range? \n",
    "\n",
    "posting_cutoff = False # Drops usernames with fewer than (mean - 0.5std) total posts (for mcmc p-val improvement)\n",
    "\n",
    "report_sample_size = False # simple reporting feature for sample size\n",
    "\n",
    "report_sm_disq = False # reports #/% of subjects disqualified for refusing to share social media data\n",
    "\n",
    "share_sm_disq_fnames = [\"../../data/final/depressioninstagram_round1.csv\",\n",
    "                       \"../../data/final/depressioninstagram_round2.csv\",\n",
    "                       \"../../data/final/controlinstagram.csv\"]\n",
    "\n",
    "path_head = '/'.join(['data-files',condition,platform])+'/'\n",
    "\n",
    "load_from = None # None, file, pickle :: loads masters from pickled file, or just masters from csv\n",
    "\n",
    "write_to = None # None, file, pickle  :: writes initial data (prepare_raw_data()) to csv or pickle\n",
    "\n",
    "final_pickle = True # pickles entire data dict after all masters are created\n",
    "\n",
    "make_hourly_plot = False # make plot of aggregated hourly posting for target vs control?\n",
    "\n",
    "run_master = True \n",
    "run_subsets = True\n",
    "run_separate_pca = True\n",
    "\n",
    "action_params = {\n",
    "    'create_master': True, \n",
    "    'save_to_file' : False, \n",
    "    'density' : False, \n",
    "    'ml' : False, \n",
    "    'nhst' : False, \n",
    "    'corr' : False, \n",
    "    'print_corrmat' : False,\n",
    "    'tall_plot': False, # (random forest) do we have enough variables to require a taller plot window?\n",
    "    'acc_avg': 'binary', # (ml) accuracy score basis for f1, precision, etc. \"binary\" or \"weighted\"\n",
    "    'use_ttest':False, # (nhst) compute t-tests for each variable?\n",
    "    'best_pca':10, # (pca) how many components to use?\n",
    "    'show_pca_comp_plot':True, # (pca) show scree plot?\n",
    "    'rf_imp_cutoff':.015, # (random forest) show feats with importance >= X\n",
    "    'rf_imp_subset':10 # (random forest) show top X important feats\n",
    "}\n",
    "\n",
    "params = define_params(condition, test_name, test_cutoff, impose_test_cutoff,\n",
    "                       platform, platform_long, fields, photos_rated, has_test, additional_data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "if not load_from: # construct data dict from raw data\n",
    "    data = make_data_dict(params, condition, test_name, conn)\n",
    "    prepare_raw_data(data, platform, params, conn, gb_types, condition, \n",
    "                     periods, turn_points, posting_cutoff,\n",
    "                     additional_data=additional_data, include_filter=include_filter, limit_date_range=limit_date_range)\n",
    "else: # load data dict from saved file\n",
    "    data = pickle.load(open(\"{path}{cond}_{pl}_data.p\".format(path=path_head,\n",
    "                                                              cond=condition,\n",
    "                                                              pl=platform), \"rb\" ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "if run_master:\n",
    "    master = data['master']\n",
    "    target = data['target']['gb']\n",
    "    control = data['control']['gb'] \n",
    "    report = 'MAIN'\n",
    "    \n",
    "    if action_params['create_master']:\n",
    "        master['model'] = {}\n",
    "\n",
    "    for gb_type in ['created_date']:\n",
    "\n",
    "        master_actions(master, target, control, condition, platform, \n",
    "                       params, gb_type, report, action_params, clfs, additional_data, posting_cutoff)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "use_pca = False # should models be fit using orthogonal pca components?\n",
    "\n",
    "if run_subsets:\n",
    "    action_params['compare_filters'] = False\n",
    "    for period in ['before']:\n",
    "        if action_params['create_master']:\n",
    "            data['master'][period] = {}\n",
    "\n",
    "        for turn_point in ['from_diag']:    \n",
    "            if action_params['create_master']:\n",
    "                data['master'][period][turn_point] = {}\n",
    "\n",
    "            master = data['master'][period][turn_point]\n",
    "            target = data['target'][period][turn_point]['gb']\n",
    "            control = data['control']['gb'] \n",
    "            report = '{}_{}'.format(period,turn_point)\n",
    "\n",
    "            if action_params['create_master']:\n",
    "                master['model'] = {}\n",
    "\n",
    "            for gb_type in ['created_date']:\n",
    "                print 'Reporting for: SUBSETS'\n",
    "                print 'Period: {}  Focus: {}  Groupby: {}'.format(period.upper(), turn_point.upper(), gb_type.upper())\n",
    "                # merge target, control, into master\n",
    "                master_actions(master, target, control, condition,\n",
    "                               platform, params, gb_type, report,\n",
    "                               action_params, clfs, additional_data, posting_cutoff,\n",
    "                               use_pca=use_pca) # using PCA?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data['target']['all']['has_face'] = None\n",
    "data['target']['all']['face_ct'] = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data['control']['all']['has_face'] = None\n",
    "data['control']['all']['face_ct'] = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print 'target shape:', data['target']['all'].shape\n",
    "print 'control shape:', data['control']['all'].shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# imports for face detect\n",
    "import cv2\n",
    "import os\n",
    "import sys\n",
    "from string import Template\n",
    "import urllib\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "from os import listdir\n",
    "from os.path import isfile, join\n",
    "import filecmp\n",
    "from IPython.display import Image, display\n",
    "import pickle\n",
    "import statsmodels.api as sm\n",
    "from sklearn.metrics import mean_squared_error"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def find_faces(row, pop, condition):\n",
    "    url = row['url']\n",
    "    path_url = url.split('//scontent.cdninstagram.com/')[1].replace('/','___') + \".jpg\"\n",
    "    savepath = '/'.join([\"photos\",condition,pop,\"{}\".format(path_url)])\n",
    "    \n",
    "    if not isfile(savepath):\n",
    "        urllib.urlretrieve(url, savepath)\n",
    "    \n",
    "    row['has_face'], row['face_ct'] = face_detect(savepath)\n",
    "\n",
    "    return row"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def face_detect(img):\n",
    "    ''' Detects faces in a photograph using Haar cascades '''\n",
    "    ''' Modified from source: https://gist.github.com/dannguyen/cfa2fb49b28c82a1068f '''\n",
    "    \n",
    "    # first argument is the haarcascades path\n",
    "    face_cascade_path = \"../../haarcascade_frontalface_default.xml\"\n",
    "    face_cascade = cv2.CascadeClassifier(face_cascade_path)\n",
    "\n",
    "    # profiles didn't work better than default (and also didn't catch things default missed)\n",
    "    #profile_cascade_path = \"../../haarcascade_profileface.xml\"\n",
    "    #profile_cascade = cv2.CascadeClassifier(profile_cascade_path)\n",
    "\n",
    "    scale_factors = [1.05, 1.4] # you played around and found these two runs to be the best\n",
    "    min_neighbors = 4 # you experimented between 1 and 5 here, 4 seemed best.\n",
    "    min_size = (20,20) # don't go higher than this, many tiny faces\n",
    "    flags = cv2.CASCADE_SCALE_IMAGE\n",
    "\n",
    "    image = cv2.imread(img)\n",
    "\n",
    "    found_face = False\n",
    "    face_ct = 0\n",
    "    \n",
    "    for scale in scale_factors:\n",
    "        if not found_face:\n",
    "            scale_factor = scale\n",
    "            for cascade, view in [(face_cascade,'straight')]:#[(face_cascade,'straight'),(profile_cascade,'profile')]:\n",
    "                faces = cascade.detectMultiScale(image, \n",
    "                                                scaleFactor = scale_factor, \n",
    "                                                minNeighbors = min_neighbors,\n",
    "                                                minSize = min_size, \n",
    "                                                flags = flags)\n",
    "                if len(faces) != 0:\n",
    "                    found_face = True\n",
    "                    face_ct += len(faces)\n",
    "\n",
    "                    for ( x, y, w, h ) in faces:\n",
    "                        cv2.rectangle(image, (x, y), (x + w, y + h), (255, 255, 0), 2)\n",
    "                        savedir = '/'.join(['photos',condition,pop,'detected',view])\n",
    "                        savefname = img.split('/')[-1] + \".jpg\"\n",
    "                        savepath = join(savedir,savefname)\n",
    "                        cv2.imwrite(savepath, image)\n",
    "                    \n",
    "    return(found_face, face_ct)            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def show_photo(url, pop, pred_face, doPrint=True):\n",
    "    ''' Displays photo in Jupyter notebook for face assessment '''\n",
    "    \n",
    "    if pred_face:\n",
    "        path_head = \"photos/depression/{}/detected/straight\".format(pop)\n",
    "    else:\n",
    "        path_head = \"photos/depression/{}\".format(pop)\n",
    "        \n",
    "    new_fname = url.split('//scontent.cdninstagram.com/')[1].replace('/','___') + \".jpg\"\n",
    "    \n",
    "    if doPrint:\n",
    "        print new_fname\n",
    "    \n",
    "    path = join(path_head, new_fname)\n",
    "    \n",
    "    try:\n",
    "        display(Image(filename=path))\n",
    "    except:\n",
    "        display(Image(filename=path+'.jpg'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running face detection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## this is the main function for pulling down photos from url and finding faces in them!\n",
    "## it takes a long time to run all photos in the depression set (hours)\n",
    "\n",
    "for pop in ['target','control']:\n",
    "    #pop = 'control'\n",
    "    data[pop]['all'] = data[pop]['all'].apply(find_faces, axis=1, args=(pop, condition))\n",
    "\n",
    "pickle.dump( data, open( \"with_face_data_{}.p\".format(condition), \"wb\" ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## updates db with has_face and face_ct fields\n",
    "## (only do this once)\n",
    "\n",
    "#with conn:\n",
    "#    cur = conn.cursor()\n",
    "#    cur.execute('alter table meta_ig add column has_face int')\n",
    "#    cur.execute('alter table meta_ig add column face_ct int')\n",
    "#    conn.commit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "## uploads face data into meta_ig\n",
    "\n",
    "tmp = data['target']['all'][['has_face','face_ct','url']].copy()\n",
    "tmp.has_face = tmp.has_face.astype(int)\n",
    "tmp.url = tmp.url.astype(str)\n",
    "\n",
    "tmp2 = data['control']['all'][['has_face','face_ct','url']].copy()\n",
    "tmp2.has_face = tmp2.has_face.astype(int)\n",
    "tmp2.url = tmp2.url.astype(str)\n",
    "\n",
    "tmp3 = pd.concat([tmp,tmp2])\n",
    "\n",
    "tups = [tuple(x) for x in tmp3[['has_face','face_ct','url']].values]\n",
    "\n",
    "q = 'update meta_ig set has_face = ?, face_ct = ? where url = ?'\n",
    "with conn:\n",
    "    cur = conn.cursor()\n",
    "    cur.executemany(q, tups)\n",
    "    conn.commit() "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Creating validation sample sets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## get all data with face ratings for condition\n",
    "\n",
    "condition = 'ptsd'\n",
    "q = 'select username, url, has_face, face_ct, d_from_diag_{cond} as diagdate, d_from_susp_{cond} as suspdate from meta_ig where has_face is not null'.format(cond=condition)\n",
    "a = pd.read_sql_query(q, conn)\n",
    "a = a.ix[a.url.isin(tmp3.url),:]\n",
    "print 'Face ratings for {}'.format(condition), a.shape\n",
    "\n",
    "## break down into control/target claim-yes-face/claim-no-face\n",
    "\n",
    "cmask = a.diagdate.isnull()\n",
    "tmask = a.diagdate.notnull()\n",
    "facemask = a.has_face\n",
    "\n",
    "cnotf = a.ix[cmask & ~facemask,:].copy()\n",
    "# control face predicted\n",
    "chasf = a.ix[cmask & facemask,:].copy()\n",
    "\n",
    "# target no face predicted\n",
    "tnotf = a.ix[tmask & ~facemask,:].copy()\n",
    "# target face predicted\n",
    "thasf = a.ix[tmask & facemask,:].copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print \"control pop:\", cnotf.shape[0] + chasf.shape[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print \"target pop:\", tnotf.shape[0] + thasf.shape[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "## get 100-observation sample from one at a time\n",
    "## set current to determine which group you're reviewing\n",
    "\n",
    "current = thasf\n",
    "\n",
    "current['actual_has_face'] = None\n",
    "current['actual_face_ct'] = None\n",
    "current['actual_face_obstruct'] = None\n",
    "current['actual_face_toosmall'] = None\n",
    "current_sample = current.sample(100).reset_index(drop=True)\n",
    "\n",
    "samp = current_sample"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Assessing faces manually"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## this block actually displays a photo to assess\n",
    "\n",
    "# set pop based on target or control\n",
    "pop = 'target'\n",
    "# set pred_face based on whether you're looking at a sample that claimed faces vs claimed no faces\n",
    "pred_face = True\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "# go one at a time!\n",
    "url_num = 99\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "show_photo(samp.url[url_num], pop, pred_face)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "## report results from photo\n",
    "\n",
    "actual_has_face      = 0\n",
    "actual_face_count    = 0\n",
    "\n",
    "actual_face_obstruct = 0 # these are counts\n",
    "actual_face_toosmall = 0 # these are counts\n",
    "\n",
    "\n",
    "samp.ix[url_num, 'actual_has_face'] = bool(actual_has_face)\n",
    "if not bool(actual_has_face):\n",
    "    samp.ix[url_num, 'actual_face_ct'] = 0\n",
    "else:\n",
    "    samp.ix[url_num, 'actual_face_ct'] = actual_face_count\n",
    "    samp.ix[url_num, 'actual_face_obstruct'] = actual_face_obstruct\n",
    "    samp.ix[url_num, 'actual_face_toosmall'] = actual_face_toosmall"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bringing together validation sets for analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# brings all four sample classes together [ target/control, pred_face/pred_no_face ]\n",
    "allf = pd.concat([thasf, tnotf])\n",
    "allf = pd.concat([allf, chasf])\n",
    "allf = pd.concat([allf, cnotf])\n",
    "\n",
    "allf.shape\n",
    "\n",
    "aa = a.ix[a.url.isin(allf.url),:].reset_index(drop=True).copy()\n",
    "allf2 = pd.merge(allf,aa, on=['url','has_face','face_ct','diagdate','username'], how='inner')\n",
    "\n",
    "# just a check, shouldn't change \n",
    "print 'before dropna:', allf2.shape\n",
    "print 'after dropna:', allf2.dropna(subset=['username']).shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def get_face_stats(subset, when='before', gb_type='created_date'):\n",
    "    ''' Reports descriptive stats for face detection in given population subset '''\n",
    "    \n",
    "    if subset == 'main':\n",
    "        aset = data['master'][gb_type]\n",
    "        subdf = 'main'\n",
    "    else:\n",
    "        subdf = \"from_{}\".format(subset)\n",
    "        col = \"{}_{}\".format(when, subset)\n",
    "        turn = '{}date'.format(subset)\n",
    "\n",
    "        aset = data['master'][when][subdf][gb_type].reset_index()\n",
    "\n",
    "        \n",
    "    # masks for all depression (target + control)\n",
    "    t_mask = aset.target==1\n",
    "    c_mask = aset.target==0\n",
    "    \n",
    "    print 'face ct avg (target):', aset.ix[t_mask,'face_ct'].mean()\n",
    "    print 'face ct std (target):', aset.ix[t_mask,'face_ct'].std()\n",
    "    print\n",
    "    print 'face ct avg (control):', aset.ix[c_mask,'face_ct'].mean()\n",
    "    print 'face ct std (control):', aset.ix[c_mask,'face_ct'].std()\n",
    "    \n",
    "    \n",
    "    #aset = a.ix[a.username.isin(bdate.username) & (a[turn].isnull() | (a[turn]<0)),:]\n",
    "\n",
    "    print\n",
    "    print 'For all data (not just validation samples...):'\n",
    "    print '{} set size: {}'.format(subdf,aset.shape[0])\n",
    "    print 'target set size:', aset.ix[t_mask,:].shape[0]\n",
    "    print 'control set size:', aset.ix[c_mask,:].shape[0]\n",
    "    print\n",
    "    print\n",
    "    \n",
    "    print 'Prop HAS FACE for TARGET:', round(aset.ix[t_mask, 'has_face'].mean(), 3)\n",
    "    print 'Std prop HAS FACE for TARGET:', round(aset.ix[t_mask, 'has_face'].std(), 3)\n",
    "    print\n",
    "    print 'Prop HAS FACE for CONTROL:', round(aset.ix[c_mask, 'has_face'].mean(), 3)\n",
    "    print 'Std prop HAS FACE for CONTROL:', round(aset.ix[c_mask, 'has_face'].std(), 3)\n",
    "\n",
    "    print\n",
    "    print\n",
    "\n",
    "    print 'Considering all photos with at least one face...'\n",
    "    print\n",
    "    print 'Mean FACE CT for TARGET:', round(aset.ix[t_mask & (aset.face_ct>0), 'face_ct'].mean(), 3)\n",
    "    print 'STD FACE CT for TARGET:', round(aset.ix[t_mask & (aset.face_ct>0), 'face_ct'].std(), 3)\n",
    "    print\n",
    "    print 'Mean FACE CT for CONTROL:', round(aset.ix[c_mask & (aset.face_ct>0), 'face_ct'].mean(), 3)\n",
    "    print 'STD FACE CT for CONTROL:', round(aset.ix[c_mask & (aset.face_ct>0), 'face_ct'].std(), 3)\n",
    "\n",
    "    print\n",
    "    print\n",
    "\n",
    "    print 'ttest for has_face:'\n",
    "    tout = ttest(aset.ix[t_mask, 'has_face'], aset.ix[c_mask, 'has_face'])\n",
    "    print 't = {}, p = {}'.format(tout.statistic, tout.pvalue)\n",
    "    print\n",
    "\n",
    "    print 'ttest for face_ct:'\n",
    "    tout = ttest(aset.ix[t_mask & (aset.face_ct>0), 'face_ct'], aset.ix[c_mask & (aset.face_ct>0), 'face_ct'])\n",
    "    print 't = {}, p = {}'.format(tout.statistic, tout.pvalue)\n",
    "    \n",
    "    print 'ttest for face_ct all photos:'\n",
    "    tout = ttest(aset.ix[t_mask, 'face_ct'], aset.ix[c_mask, 'face_ct'])\n",
    "    print 't = {}, p = {}'.format(tout.statistic, tout.pvalue)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "subset = \"diag\" # takes values: from_diag, from_susp, main\n",
    "\n",
    "get_face_stats(subset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# masks for samples\n",
    "mask_t = allf2.diagdate.notnull()\n",
    "mask_c = ~mask_t"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## checking to make sure we don't have any high-face-posting outliers in depression group that would explain diff\n",
    "\n",
    "(allf2.ix[mask_t, :]\n",
    "      .groupby('username')\n",
    "      .agg({'face_ct':'mean'})['face_ct']\n",
    "      .plot(kind='hist', \n",
    "            bins=10, \n",
    "            xlim=(0,15), \n",
    "            ylim=(0,15),\n",
    "            title=\"Target group: avg face count per user\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## actually, we have a high-posting outlier in control!  without it results would be same but more extreme\n",
    "\n",
    "(allf2.ix[mask_c, :]\n",
    "      .groupby('username')\n",
    "      .agg({'face_ct':'mean'})['face_ct']\n",
    "      .plot(kind='hist', \n",
    "            bins=50, \n",
    "            xlim=(0,15), \n",
    "            ylim=(0,25),\n",
    "            title=\"Control group: avg face count per user\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print 'For validation samples:'\n",
    "print\n",
    "\n",
    "print 'TARGET :: true negative :: has_face:', 1 - allf2.ix[mask_t & ~allf2.has_face, 'actual_has_face'].mean()\n",
    "print 'CONTROL :: true negative :: has_face:', 1 - allf2.ix[mask_c & ~allf2.has_face, 'actual_has_face'].mean()\n",
    "\n",
    "print \n",
    "\n",
    "print 'TARGET :: true positive :: has_face:', allf2.ix[mask_t & (allf2.has_face), 'actual_has_face'].mean()\n",
    "print 'CONTROL :: true positive :: has_face:', allf2.ix[mask_c & allf2.has_face, 'actual_has_face'].mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## Looking at face ct accuracy differences between target and controls\n",
    "\n",
    "t_ct = allf2.ix[mask_t , 'face_ct'].values\n",
    "t_act_ct = allf2.ix[mask_t , 'actual_face_ct'].values\n",
    "\n",
    "print 'TARGET face ct vs actual face ct MEAN DIFF ct-act_ct', np.mean( t_ct - t_act_ct)\n",
    "print 'TARGET face ct vs actual face ct STD', np.std(t_ct - t_act_ct)\n",
    "mse = mean_squared_error(t_ct, t_act_ct)\n",
    "print 'TARGET face ct vs actual face ct MSE', mse\n",
    "\n",
    "print\n",
    "print\n",
    "\n",
    "c_ct = allf2.ix[mask_c , 'face_ct'].values\n",
    "c_act_ct = allf2.ix[mask_c , 'actual_face_ct'].values\n",
    "\n",
    "print 'CONTROL face ct vs actual face ct MEAN DIFF ct-act_ct', np.mean( c_ct - c_act_ct)\n",
    "print 'CONTROL face ct vs actual face ct STD', np.std(c_ct - c_act_ct)\n",
    "mse = mean_squared_error(c_ct, c_act_ct)\n",
    "print 'CONTROL face ct vs actual face ct MSE', mse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## logistic regression with has_face, face_ct as predictors\n",
    "\n",
    "a['target'] = None\n",
    "a.ix[a.diagdate.notnull(), 'target'] = 1\n",
    "a.ix[a.diagdate.isnull(), 'target'] = 0\n",
    "\n",
    "dm = pd.DataFrame(a[['has_face','face_ct']].copy())\n",
    "dm = sm.add_constant(dm)\n",
    "logit = sm.Logit(a.target.astype(int).values, dm)\n",
    "# fit the model\n",
    "result = logit.fit()\n",
    "print result.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "## Makeup samples for controls with missing images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## use this block for samples where no face was found, and where some images just fail to load\n",
    "\n",
    "#new_samp_size = 100-total.shape[0] \n",
    "#sample2 = [insert_current_df_like_cnotf...].sample(new_samp_size).reset_index(drop=True)\n",
    "\n",
    "\n",
    "\n",
    "## in these cases, you took a new sample to make up the no-load photos, so each sample total still summed to 100\n",
    "## this code block concats the earlier sample and the makeup sample set.\n",
    "\n",
    "# dhasf_sample = samp\n",
    "# dhasf_sample['pop'] = 'target'\n",
    "\n",
    "#tmp = dnotf_sample.copy()\n",
    "#tmp.dropna(subset=['actual_has_face'], inplace=True)\n",
    "#tmp = pd.concat([tmp, samp], axis=0)\n",
    "#tmp.dropna(subset=['actual_has_face']).shape\n",
    "\n",
    "## now pickle file to disk\n",
    "\n",
    "#pickle.dump(allf,open('allf.p','wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# when you originally did this for depression, allf2 ended up being a corrected version of allf validation sample data.  \n",
    "# 'allf2_USE_THIS_ALLF.p' is the right file to use if you need the validation samples from depression.\n",
    "\n",
    "#pickle.dump(allf2,open('allf2_USE_THIS_ALLF.p','wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}