In [70]:
from glob import glob
import seaborn as sns
import nibabel as nib
import numpy as np 
from nilearn import image, plotting
import json
import os.path as osp
import pandas as pd
from matplotlib import pyplot as plt
import statsmodels
import statsmodels.formula.api as smf
%matplotlib inline
sns.set(rc={'axes.facecolor':'lightgray', 'figure.facecolor':'white'})

First with TBSS images


In [258]:
mds = glob('/home/grg/tbss_subjects/*L2.nii.gz')
skel = np.array(nib.load('/home/grg/tbss_fmodel/stats/mean_FA_skeleton.nii.gz').dataobj)
ages = json.load(open('/home/grg/spm/data/age.json'))
apoe = json.load(open('/home/grg/spm/data/apoe_groups.json'))
gender = json.load(open('/home/grg/spm/data/genders.json'))
educ = json.load(open('/home/grg/spm/data/educyears.json'))

data = []
for each in mds:    
    md = np.array(nib.load(each).dataobj)
    mean = md[skel>0.01].mean()
    subject = osp.basename(each).split('_')[1]
    print subject
    data.append((subject, mean, ages[subject[:5]]/365.25, apoe[subject], educ[subject], gender[subject]))


10117
13235
13214
10737
12502
10900
66164
12841
10166
66048
11196
10703
12479
55539
12172
66309
10251
11184
13242
10298
10613
11620
12133
11592
21039
11007
10419
12767
10811
12445
10988
13345
10029
12484
10282
10822
11222
12711
10015
12975
10035
10024
10576
12765
44046
21051
12079
11691
10226
10170
12947
11630
10436
11038
10144
11882
12245
21002
10855
11407
12636
10070
10156
11063
10036
11646
11133
11436
21130
10846
12483
10065
12582
11144
10118
11225
13083
12067
10541
77179
10418
12699
10013
12304
66085
11245
10042
10180
12271
10496
13070
11254
10322
12930
10538
13035
77252
11939
10645
11137
10334
55216
11195
10034
11291
13061
12391
11262
55630
10212
11737
12787
11850
13169
10099
11019
11045
44660
11550
10225
10326
13138
12279
12963
10330
44057
12874
11414
12812
10367
11152
66125
77195
66492
13215
21092
10053
10239c
---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-258-94f9ab392661> in <module>()
     14     subject = osp.basename(each).split('_')[1]
     15     print subject
---> 16     data.append((subject, mean, ages[subject[:5]]/365.25, apoe[subject], educ[subject], gender[subject]))

KeyError: '10239c'

In [274]:
df


Out[274]:
roi age apo educyears gender group
subject
10117 0.000988 63.008898 2 17 0 0
13235 0.000961 56.177960 0 8 0 0
13214 0.000964 48.290212 2 13 0 0
10737 0.001007 49.180014 0 15 0 0
12502 0.001055 61.765914 4 12 1 2
10900 0.001102 65.711157 2 10 0 0
66164 0.000960 65.185489 1 18 1 1
12841 0.000877 48.418891 0 8 0 0
10166 0.000991 48.947296 2 11 0 0
66048 0.001006 68.123203 2 18 1 0
11196 0.000952 56.032854 3 18 1 1
10703 0.000964 48.465435 0 18 1 0
12479 0.000990 48.197125 0 10 0 0
55539 0.001031 59.008898 1 8 0 1
12172 0.000955 58.124572 4 8 0 2
66309 0.000990 61.982204 4 12 0 2
10251 0.000919 54.617385 2 15 0 0
11184 0.001125 64.372348 2 18 1 0
13242 0.000971 56.971937 2 9 0 0
10298 0.000923 48.714579 4 18 0 2
10613 0.000989 51.403149 4 17 0 2
11620 0.000932 51.271732 2 18 1 0
12133 0.000918 51.805613 4 11 0 2
11592 0.000954 55.770021 4 12 1 2
21039 0.000957 49.292266 4 18 1 2
11007 0.000980 57.527721 3 17 1 1
10419 0.000915 50.904860 0 8 0 0
12767 0.001025 49.804244 2 10 0 0
10811 0.000980 57.475702 3 18 0 1
12445 0.000913 58.162902 0 17 0 0
... ... ... ... ... ... ...
11262 0.000963 64.659822 0 8 1 0
55630 0.000944 54.360027 4 14 1 2
10212 0.000959 64.933607 4 17 1 2
11737 0.000994 65.349760 0 12 1 0
12787 0.000943 65.478439 0 10 0 0
11850 0.000977 53.242984 3 15 1 1
13169 0.000931 63.939767 2 8 0 0
10099 0.001035 51.663244 3 11 0 1
11019 0.000970 48.657084 3 17 0 1
11045 0.000972 64.851472 0 17 1 0
44660 0.000970 57.902806 4 10 0 2
11550 0.000989 53.037645 0 18 1 0
10225 0.001026 61.921971 3 10 0 1
10326 0.001038 50.869268 2 18 1 0
13138 0.001016 49.587953 0 11 0 0
12279 0.000969 55.488022 2 12 1 0
12963 0.000997 47.707050 4 18 1 2
10330 0.000878 49.779603 2 15 1 0
44057 0.001043 66.844627 3 10 1 1
12874 0.001012 51.737166 3 10 0 1
11414 0.000975 53.108830 0 17 0 0
12812 0.000954 45.681040 3 8 0 1
10367 0.000976 53.946612 4 8 1 2
11152 0.000938 61.686516 2 8 0 0
66125 0.000992 52.027379 0 18 0 0
77195 0.000964 60.271047 0 12 0 0
66492 0.000913 52.391513 1 11 0 1
13215 0.000959 47.411362 2 11 0 0
21092 0.000981 52.147844 0 17 1 0
10053 0.000984 50.116359 2 18 0 0

137 rows × 6 columns


In [283]:
def load_md(fp):
    df = pd.read_excel(fp)
    subjects = df['subject'].tolist()
    df = df.set_index('subject')

    #d = {'apoe23':0, 'apoe24':1, 'apoe33':0, 'apoe34':1, 'apoe44':2}
    d = {0:0, 1:1, 2:0, 3:1, 4:2}

    ht = {0:1, 1:2, 2:1, 3:2, 4:1}
    groups = [d[e] for e in df['apo'].values]

    df['group'] = pd.Series(groups, index=subjects)
    return df

In [313]:
df_dartel = load_md('/tmp/mean_MD_dartel.xls').dropna()
df_tbss = load_md('/tmp/mean_MD.xls').dropna()

%run /home/grg/git/alfa/nilearn-helper.py
%run /home/grg/git/alfa/roicollect.py
#sns.boxplot(x='groups', y='mean_MD', data=df, showfliers=False)

#groups = get_groups(df, )
df_dartel['roi'] = correct(df_dartel)
df_tbss['roi'] = correct(df_tbss)
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
set_figaxes(df_tbss, ylim=[df_dartel['roi'].min(), df_tbss['roi'].max()])

plot_region(df_dartel, '', ['All'], order=1, ax=ax)
plot_region(df_tbss, '', ['All'], order=1, ax=ax, c='toto')
#sns.boxplot(x='group', y='roi', data=df, showfliers=False)#, ax=ax)
#df = df.dropna()
#boxplot_region(df, ['HO','HT','NC'])


Model used for correction: roi ~ 1 + gender + educyears + apo
Model used for correction: roi ~ 1 + gender + educyears + apo
Region: 
Region: 

In [301]:



dartel quadratic aic= -8995.86168923 bic= -8983.20588888
dartel linear aic= -8996.34771512 bic= -8987.91051489
tbss quadratic aic= -2339.36106155 bic= -2330.60111877
tbss linear aic= -2340.32947673 bic= -2334.48951488

Now with DARTEL images


In [264]:
mds = sorted(glob('/home/grg/dartel_csf.5/rswr*_MD_t1space_wo_csf_nohdr.nii'))
#skel = sorted(glob('/home/grg/data/ALFA_DWI/*/T1/*_mabonlm_nobias_spm_c2.nii'))
ages = json.load(open('/home/grg/spm/data/age.json'))
apoe = json.load(open('/home/grg/spm/data/apoe_groups.json'))
educ = json.load(open('/home/grg/spm/data/educyears.json'))
gender = json.load(open('/home/grg/spm/data/genders.json'))

data = []
for each in mds:    
    #md = np.array(nib.load(each).dataobj)
    subject = osp.basename(each).split('_')[0][4:]
    print subject
    skelfp = '/home/grg/data/templates/MNI_AAL/MNI_T1_seg_2.nii.gz'
    # glob('/home/grg/data/ALFA_DWI/%s*/T1/%s*_mabonlm_nobias_spm_c2.nii'%(subject, subject))[0]
    #image.resample_to_img(each, skelfp).to_filename('/tmp/teeesst.nii.gz')
    md = np.array(nib.load(each).dataobj)
    skel = np.array(nib.load(skelfp).dataobj)
    mean = md[skel>0.01].mean()
    print mean
    data.append((subject, mean, ages[subject[:5]]/365.25, apoe[subject], educ[subject], gender[subject]))

df = pd.DataFrame(data, columns=['subject', 'roi', 'age', 'apo', 'educyears', 'gender'])
df.to_excel('/tmp/mean_MD_dartel.xls')
df.head()


10013
0.000803078
10015
0.000903842
10016
0.000815249
10023
0.000846521
10024
0.000880902
10025
0.000821792
10026
0.000839239
10028
0.000862139
10029
0.000847982
10034
0.000842253
10035
0.000825732
10036
0.000886464
10038
0.000871669
10040
0.000909516
10041
0.000815894
10042
0.000862737
10049
0.000889316
10051
0.000886537
10052
0.000859474
10053
0.000858774
10056
0.000868197
10065
0.000847915
10070
0.000861453
10071
0.000805527
10081
0.000880744
10090
0.000858479
10096
0.00094524
10099
0.000885283
10102
0.000855912
10106
0.000886683
10108
0.000848911
10117
0.000857764
10118
0.000866308
10134
0.000836692
10144
0.000868438
10150
0.000836771
10151
0.000845446
10156
0.000833625
10158
0.000832036
10160
0.00086361
10162
0.000836652
10166
0.000880331
10170
0.000830705
10178
0.000930228
10180
0.000820001
10182
0.000849002
10199
0.00085866
10200
0.000876952
10212
0.000868673
10213
0.000832982
10217
0.000815728
10225
0.000822936
10226
0.000893113
10235
0.000871713
10242
0.000885162
10245
0.000840485
10248
0.000875212
10251
0.0008375
10253
0.000827656
10259
0.000891471
10263
0.00083446
10265
0.000895105
10282
0.000818119
10298
0.000822091
10308
0.000841673
10313
0.000926082
10317
0.000887552
10319
0.000878239
10322
0.000873462
10324
0.000852455
10325
0.000927532
10326
0.000854905
10329
0.000871277
10330
0.000816159
10334
0.000823209
10338
0.000883884
10346
0.000837007
10354
0.000863423
10361
0.000858203
10362
0.000871598
10365
0.000857514
10367
0.000857538
10370
0.000883036
10385
0.000814445
10393
0.000867216
10396
0.000788737
10397
nan
10416
0.000849327
10417
0.000837536
10418
0.000831172
10419
0.000863203
10426
0.000851669
10433
0.000841835
10436
0.000840087
10450
0.000815097
10453
0.000809564
10461
0.000860912
10463
0.00083014
10482
0.000849566
10486
0.000858417
10493
0.00083551
10496
0.000814605
10504
0.000881571
10515
0.000836754
10518
0.00081406
10522
0.000870492
10528
0.00085834
10530
0.000824405
10536
0.000842629
10538
0.00086406
10541
0.000887001
10550
0.000870534
10551
nan
10563
0.000878248
10576
0.00088693
10577
0.000925958
10593
0.000910721
10613
0.000856988
10630
0.000824066
10634
0.000823372
10645
0.000893774
10657
0.000870429
10668
0.000839631
10678
0.000892289
10682
0.000819231
10692
0.000857848
10693
0.000885591
10696
0.000857484
10697
0.000841914
10701
0.000759166
10703
0.000891
10724
0.00091561
10725
0.000839691
10735
0.0008809
10737
0.00087965
10741
0.000888054
10744
0.000837911
10750
0.000840471
10756
0.000837227
10778
0.000880887
10787
0.000866845
10794
0.000897395
10809
0.000855505
10811
0.000820348
10821
0.000807015
10822
0.000854323
10841
0.000920931
10846
0.000848818
10850
0.000882899
10855
0.000898095
10858
0.000818417
10870
0.000834803
10881
0.000827551
10894
0.000831789
10900
0.000886801
10901
0.000836145
10942
0.000868793
10944
0.000849423
10946
0.00087548
10972
0.000873861
10988
0.000865548
11007
0.000855979
11012
0.000870616
11015
0.000852212
11019
0.000870582
11030
0.000870829
11038
0.000846269
11042
0.00088895
11045
0.000867761
11047
0.000879143
11048
0.000857285
11054
0.000874427
11063
0.000854922
11092
0.000844614
11114
0.000859188
11127
0.000849084
11133
0.000865082
11136
0.000830626
11137
0.000919367
11139
0.000881791
11144
0.000852473
11152
0.000840555
11156
0.000864925
11180
0.000822953
11184
0.000880764
11191
0.000836633
11195
0.00079803
11196
0.000842491
11201
0.000876638
11213
0.000860432
11219
0.000831721
11222
0.000901756
11225
0.000822491
11245
0.000837986
11247
0.000848004
11252
0.000847435
11254
0.000894801
11262
0.000852754
11264
0.000827383
11291
0.000925339
11292
0.000804618
11305
0.000840983
11323
0.000844022
11327
0.00085006
11351
0.00082688
11355
0.000880184
11360
0.000823172
11383
nan
11387
0.000919015
11407
0.000868036
11414
0.000884048
11415
0.000818633
11416
0.000875415
11426
0.000912242
11436
0.000910451
11458
0.000873161
11461
0.000869482
11474
0.000838487
11478
0.000903117
11514
0.000846092
11540
0.0008127
11550
0.000865559
11552
0.000932593
11561
0.000893035
11583
0.000846097
11590
0.000858796
11592
0.00087449
11593
0.000861059
11597
0.000862032
11610
0.00087905
11614
0.000862041
11620
0.000876048
11630
0.000825442
11638
0.000870363
11641
0.00081197
11646
0.000881685
11656
0.000857856
11658
0.000844516
11679
0.000805465
11686
0.000820502
11687
0.000876073
11691
0.000875805
11711
0.00089917
11721
0.000834848
11737
0.000866374
11747
0.000897756
11768
0.000905758
11796
0.00085764
11798
0.000801443
11803
0.000944602
11829
0.000865213
11830
0.000855655
11850
0.000830086
11858
0.0008833
11872
0.000847638
11874
0.000842048
11882
0.000898967
11902
0.000815613
11937
0.000866091
11939
0.000835777
11941
0.000834907
11943
0.000874748
11975
0.000831745
11979
0.000842892
12015
0.00087025
12032
0.000839607
12056
0.000850186
12067
0.000842225
12079
0.000872573
12121
0.000876717
12122
0.000837779
12125
0.000806753
12133
0.000816182
12138
0.000826135
12140
0.000823801
12172
0.000855101
12174
0.000849926
12186
0.000882416
12239
0.000835577
12244
0.00085197
12245
0.000811044
12252
0.000847572
12269
0.000863699
12271
0.000900149
12279
0.000814174
12296
0.000819316
12304
0.000849898
12308
0.000821508
12323
0.000840496
12324
0.000876705
12327
0.000869825
12331
0.000891121
12356
0.000841338
12379
0.000825718
12391
0.000868824
12399
0.000825266
12409
0.000871024
12417
0.000880626
12425
0.000820914
12445
0.000855248
12479
0.000865906
12483
0.000819882
12484
0.000937246
12493
0.000805086
12502
0.000941905
12511
0.000877988
12516
0.000846354
12548
0.000829711
12582
0.000860183
12624
0.000827215
12636
0.000853873
12637
0.000883968
12659
0.000823886
12699
0.000853635
12704
0.000858972
12711
0.000851268
12715
0.000822257
12724
0.000842111
12730
0.000828989
12765
0.000821527
12767
0.00085504
12771
0.00086105
12772
0.000857839
12778
0.000840474
12783
0.000834636
12785
0.000837121
12787
0.00086019
12810
0.000821361
12812
0.000851025
12823
0.000857035
12841
0.000803364
12858
0.000897179
12861
0.000885191
12874
0.000827058
12878
0.000835101
12893
0.00086958
12904
0.000817845
12920
0.000844054
12930
0.000856346
12941
0.000865856
12947
0.000875739
12963
0.000894047
12970
0.000860104
12975
0.000841185
12976
0.000858571
12995
0.000810847
13008
0.000818331
13035
0.000811195
13043
0.000875874
13049
0.000920375
13054
0.000848802
13059
0.000811659
13061
0.000874325
13063
0.000881445
13070
0.000889034
13075
0.000867057
13083
0.000850499
13090
0.000861596
13105
0.000825598
13118
0.000818061
13127
0.000823382
13138
0.000906963
13144
0.000885011
13151
0.000829154
13169
0.000830986
13188
0.000885455
13214
0.000850013
13215
0.000857316
13217
0.000863983
13235
0.000846562
13236
0.000830428
13238
0.000893306
13242
0.000846297
13244
0.000829598
13268
0.000859707
13293
0.000840058
13306
0.000937691
13309
0.000861033
13312
0.000845897
13322
0.000895697
13345
0.000900881
13367
0.000876396
13417
0.000848189
21002
0.000855678
21011
0.000880037
21012
0.000908612
21039
0.000860348
21042
0.000829097
21051
0.000814054
21056
0.000866138
21073
0.000866953
21084
0.000904597
21092
0.000647101
21130
0.000823698
44004
0.000898839
44043
0.000890919
44046
0.000894656
44057
0.000867392
44068
0.000854546
44091
0.000908031
44094
0.000933675
44119
0.000873
44141
0.000876416
44147
0.00089698
44151
0.000904859
44205
0.000857591
44229
0.000806248
44491
0.000845022
44632
0.000862367
44660
0.000859959
44723
0.000841149
55057
0.000854487
55152
0.000864945
55166
0.000824973
55200
0.000824301
55216
0.000860004
55297
0.000890197
55323
0.000870181
55351
0.000863536
55353
0.000851195
55355
0.000841859
55370
0.000796646
55388
0.000914243
55469
0.000832505
55483
0.000890722
55488
0.000823265
55529
0.000897414
55538
0.000834123
55539
0.000882261
55630
0.000856559
55636
0.000876986
55667
0.00079806
55708
0.000893898
55734
0.00087381
55778
0.000931003
55793
0.000845248
55854
0.000860709
66017
0.000849418
66019
0.000823797
66020
0.000866027
66026
0.000851358
66030
0.000878925
66039
0.000822526
66042
0.000858248
66048
0.00088057
66050
0.000884432
66085
0.000895819
66089
0.000872458
66094
0.00087627
66125
0.000888452
66128
0.000833028
66131
0.000813384
66133
0.000836018
66141
0.000848879
66159
0.000802958
66162
0.000839315
66164
0.000824726
66169
0.000886767
66183
0.000823087
66239
0.000846238
66240
0.000877143
66257
0.000883052
66264
0.000840644
66267
0.000876419
66268
0.000809063
66270
0.000846632
66293
0.000860941
66309
0.000880464
66312
0.000838267
66335
0.000847288
66341
0.00082065
66361
0.000880708
66368
0.000852641
66492
0.000812739
66498
0.000928375
77024
0.000897375
77027
0.000825969
77034
0.000868756
77037
0.000908991
77040
0.0008672
77044
0.00089123
77047
0.000867718
77056
0.000871077
77068
0.000931692
77076
0.000866316
77093
0.000820932
77094
0.000831077
77096
0.000826625
77109
0.000815115
77117
0.000828416
77130
0.000853106
77140
0.000853731
77151
0.000885849
77152
0.000855286
77175
0.000877659
77179
0.000836687
77188
0.000867236
77191
0.000812203
77192
0.000838295
77195
0.00088907
77217
0.000850833
77241
0.000856958
77252
0.000851567
77254
0.000877297
77263
0.000889996
Out[264]:
subject roi age apo educyears gender
0 10013 0.000803 55.414100 2 11 0
1 10015 0.000904 50.707734 1 17 0
2 10016 0.000815 48.583162 3 17 0
3 10023 0.000847 61.065024 2 18 1
4 10024 0.000881 49.401780 2 18 0

In [268]:
mod = smf.ols(formula='roi ~ 1 + age + I(age**2)', data=df_dartel)
res = mod.fit()
print 'dartel quadratic', 'aic=', res.aic, 'bic=', res.bic
mod = smf.ols(formula='roi ~ 1 + age', data=df_dartel)
res = mod.fit()
print 'dartel linear', 'aic=', res.aic, 'bic=', res.bic
mod = smf.ols(formula='roi ~ 1 + age + I(age**2)', data=df_tbss)
res = mod.fit()
print 'tbss quadratic', 'aic=', res.aic, 'bic=', res.bic
mod = smf.ols(formula='roi ~ 1 + age', data=df_tbss)
res = mod.fit()
print 'tbss linear', 'aic=', res.aic, 'bic=', res.bic


Model used for correction: roi ~ 1 + gender + educyears+ apo
Out[268]:
subject
10013    0.000803
10015    0.000906
10016    0.000815
10023    0.000847
10024    0.000882
10025    0.000821
10026    0.000839
10028    0.000863
10029    0.000848
10034    0.000843
10035    0.000825
10036    0.000885
10038    0.000873
10040    0.000911
10041    0.000818
10042    0.000864
10049    0.000889
10051    0.000886
10052    0.000860
10053    0.000860
10056    0.000869
10065    0.000846
10070    0.000863
10071    0.000803
10081    0.000882
10090    0.000859
10096    0.000944
10099    0.000884
10102    0.000859
10106    0.000886
           ...   
77024    0.000896
77027    0.000824
77034    0.000868
77037    0.000910
77040    0.000866
77044    0.000890
77047    0.000867
77056    0.000869
77068    0.000934
77076    0.000868
77093    0.000823
77094    0.000834
77096    0.000827
77109    0.000817
77117    0.000827
77130    0.000852
77140    0.000854
77151    0.000885
77152    0.000856
77175    0.000878
77179    0.000834
77188    0.000867
77191    0.000813
77192    0.000837
77195    0.000891
77217    0.000851
77241    0.000859
77252    0.000851
77254    0.000876
77263    0.000889
Name: roi, dtype: float64

In [307]:
mod = smf.ols(formula='roi ~ 1 + I(age**2)', data=df_tbss)
res = mod.fit()
res.summary()


Out[307]:
OLS Regression Results
Dep. Variable: roi R-squared: 0.118
Model: OLS Adj. R-squared: 0.112
Method: Least Squares F-statistic: 18.09
Date: Thu, 29 Jun 2017 Prob (F-statistic): 3.91e-05
Time: 18:35:50 Log-Likelihood: 1172.4
No. Observations: 137 AIC: -2341.
Df Residuals: 135 BIC: -2335.
Df Model: 1
Covariance Type: nonrobust
coef std err t P>|t| [95.0% Conf. Int.]
Intercept 0.0009 1.77e-05 51.082 0.000 0.001 0.001
I(age ** 2) 2.435e-08 5.72e-09 4.253 0.000 1.3e-08 3.57e-08
Omnibus: 5.207 Durbin-Watson: 2.119
Prob(Omnibus): 0.074 Jarque-Bera (JB): 4.734
Skew: 0.436 Prob(JB): 0.0938
Kurtosis: 3.262 Cond. No. 1.36e+04

In [ ]: