We are currently solving the following issue within our repository. Here is the issue text:
--- BEGIN ISSUE ---
Misleading (wrong?) code in Oasis VBM example
As pointed out on http://neurostars.org/p/3388/, the Oasis VBM code applies the NiftiMasker twice, and thus ends up smoothing the data twice.

It does that to do screening of low variance features. I think that the right way to do this would be to use a sklearn.feature_selection.VarianceThreshold in the pipeline.


--- END ISSUE ---

Below are some code segments, each from a relevant file. One or more of these files may contain bugs.
--- BEGIN FILES ---
Path: `examples/decoding/plot_oasis_vbm.py`
Content:
```
1 """
2 Voxel-Based Morphometry on Oasis dataset
3 ========================================
4 
5 This example uses Voxel-Based Morphometry (VBM) to study the relationship
6 between aging and gray matter density.
7 
8 The data come from the `OASIS <http://www.oasis-brains.org/>`_ project.
9 If you use it, you need to agree with the data usage agreement available
10 on the website.
11 
12 It has been run through a standard VBM pipeline (using SPM8 and
13 NewSegment) to create VBM maps, which we study here.
14 
15 Predictive modeling analysis: VBM bio-markers of aging?
16 --------------------------------------------------------
17 
18 We run a standard SVM-ANOVA nilearn pipeline to predict age from the VBM
19 data. We use only 100 subjects from the OASIS dataset to limit the memory
20 usage.
21 
22 Note that for an actual predictive modeling study of aging, the study
23 should be ran on the full set of subjects. Also, parameters such as the
24 smoothing applied to the data and the number of features selected by the
25 Anova step should be set by nested cross-validation, as they impact
26 significantly the prediction score.
27 
28 Brain mapping with mass univariate
29 -----------------------------------
30 
31 SVM weights are very noisy, partly because heavy smoothing is detrimental
32 for the prediction here. A standard analysis using mass-univariate GLM
33 (here permuted to have exact correction for multiple comparisons) gives a
34 much clearer view of the important regions.
35 
36 ____
37 
38 """
39 # Authors: Elvis Dhomatob, <elvis.dohmatob@inria.fr>, Apr. 2014
40 #          Virgile Fritsch, <virgile.fritsch@inria.fr>, Apr 2014
41 #          Gael Varoquaux, Apr 2014
42 import numpy as np
43 from scipy import linalg
44 import matplotlib.pyplot as plt
45 from nilearn import datasets
46 from nilearn.input_data import NiftiMasker
47 
48 n_subjects = 100  # more subjects requires more memory
49 
50 ### Load Oasis dataset ########################################################
51 oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=n_subjects)
52 gray_matter_map_filenames = oasis_dataset.gray_matter_maps
53 age = oasis_dataset.ext_vars['age'].astype(float)
54 
55 # print basic information on the dataset
56 print('First gray-matter anatomy image (3D) is located at: %s' %
57       oasis_dataset.gray_matter_maps[0])  # 3D data
58 print('First white-matter anatomy image (3D) is located at: %s' %
59       oasis_dataset.white_matter_maps[0])  # 3D data
60 
61 ### Preprocess data ###########################################################
62 nifti_masker = NiftiMasker(
63     standardize=False,
64     smoothing_fwhm=2,
65     memory='nilearn_cache')  # cache options
66 # remove features with too low between-subject variance
67 gm_maps_masked = nifti_masker.fit_transform(gray_matter_map_filenames)
68 gm_maps_masked[:, gm_maps_masked.var(0) < 0.01] = 0.
69 # final masking
70 new_images = nifti_masker.inverse_transform(gm_maps_masked)
71 gm_maps_masked = nifti_masker.fit_transform(new_images)
72 n_samples, n_features = gm_maps_masked.shape
73 print("%d samples, %d features" % (n_subjects, n_features))
74 
75 ### Prediction with SVR #######################################################
76 print("ANOVA + SVR")
77 # Define the prediction function to be used.
78 # Here we use a Support Vector Classification, with a linear kernel
79 from sklearn.svm import SVR
80 svr = SVR(kernel='linear')
81 
82 # Dimension reduction
83 from sklearn.feature_selection import SelectKBest, f_regression
84 
85 # Here we use a classical univariate feature selection based on F-test,
86 # namely Anova.
87 feature_selection = SelectKBest(f_regression, k=2000)
88 
89 # We have our predictor (SVR), our feature selection (SelectKBest), and now,
90 # we can plug them together in a *pipeline* that performs the two operations
91 # successively:
92 from sklearn.pipeline import Pipeline
93 anova_svr = Pipeline([('anova', feature_selection), ('svr', svr)])
94 
95 ### Fit and predict
96 anova_svr.fit(gm_maps_masked, age)
97 age_pred = anova_svr.predict(gm_maps_masked)
98 
99 # Visualization
100 # Look at the SVR's discriminating weights
101 coef = svr.coef_
102 # reverse feature selection
103 coef = feature_selection.inverse_transform(coef)
104 # reverse masking
105 weight_img = nifti_masker.inverse_transform(coef)
106 
107 # Create the figure
108 from nilearn.plotting import plot_stat_map, show
109 bg_filename = gray_matter_map_filenames[0]
110 z_slice = 0
111 from nilearn.image.resampling import coord_transform
112 affine = weight_img.get_affine()
113 _, _, k_slice = coord_transform(0, 0, z_slice,
114                                 linalg.inv(affine))
115 k_slice = np.round(k_slice)
116 
117 fig = plt.figure(figsize=(5.5, 7.5), facecolor='k')
118 weight_slice_data = weight_img.get_data()[..., k_slice, 0]
119 vmax = max(-np.min(weight_slice_data), np.max(weight_slice_data)) * 0.5
120 display = plot_stat_map(weight_img, bg_img=bg_filename,
121                         display_mode='z', cut_coords=[z_slice],
122                         figure=fig, vmax=vmax)
123 display.title('SVM weights', y=1.2)
124 
125 # Measure accuracy with cross validation
126 from sklearn.cross_validation import cross_val_score
127 cv_scores = cross_val_score(anova_svr, gm_maps_masked, age)
128 
129 # Return the corresponding mean prediction accuracy
130 prediction_accuracy = np.mean(cv_scores)
131 print("=== ANOVA ===")
132 print("Prediction accuracy: %f" % prediction_accuracy)
133 print("")
134 
135 ### Inference with massively univariate model #################################
136 print("Massively univariate model")
137 
138 # Statistical inference
139 from nilearn.mass_univariate import permuted_ols
140 neg_log_pvals, t_scores_original_data, _ = permuted_ols(
141     age, gm_maps_masked,  # + intercept as a covariate by default
142     n_perm=1000,  # 1,000 in the interest of time; 10000 would be better
143     n_jobs=1)  # can be changed to use more CPUs
144 signed_neg_log_pvals = neg_log_pvals * np.sign(t_scores_original_data)
145 signed_neg_log_pvals_unmasked = nifti_masker.inverse_transform(
146     signed_neg_log_pvals)
147 
148 # Show results
149 threshold = -np.log10(0.1)  # 10% corrected
150 
151 fig = plt.figure(figsize=(5.5, 7.5), facecolor='k')
152 
153 display = plot_stat_map(signed_neg_log_pvals_unmasked, bg_img=bg_filename,
154                         threshold=threshold, cmap=plt.cm.RdBu_r,
155                         display_mode='z', cut_coords=[z_slice],
156                         figure=fig)
157 title = ('Negative $\log_{10}$ p-values'
158          '\n(Non-parametric + max-type correction)')
159 display.title(title, y=1.2)
160 
161 signed_neg_log_pvals_slice_data = \
162     signed_neg_log_pvals_unmasked.get_data()[..., k_slice, 0]
163 n_detections = (np.abs(signed_neg_log_pvals_slice_data) > threshold).sum()
164 print('\n%d detections' % n_detections)
165 
166 show()
167
```
--- END FILES ---

Please first localize the bug based on the issue statement, and then generate a patch according to the `git diff` format fenced by three backticks.

Here is an example:

```diff
diff --git a/examples/server_async.py b/examples/server_async.py
--- a/examples/server_async.py
+++ b/examples/server_async.py
@@ -313,4 +313,4 @@
 
 
 if __name__ == "__main__":
-    asyncio.run(run_async_server("."), debug=True)
+    asyncio.run(run_async_server(), debug=True)
diff --git a/examples/server_sync.py b/examples/server_sync.py
--- a/examples/server_sync.py
+++ b/examples/server_sync.py
@@ -313,5 +313,5 @@
 
 
 if __name__ == "__main__":
-    server = run_sync_server(".")
+    server = run_sync_server()
     server.shutdown()

```

diff --git a/examples/decoding/plot_oasis_vbm.py b/examples/decoding/plot_oasis_vbm.py
--- a/examples/decoding/plot_oasis_vbm.py
+++ b/examples/decoding/plot_oasis_vbm.py
@@ -63,12 +63,7 @@
     standardize=False,
     smoothing_fwhm=2,
     memory='nilearn_cache')  # cache options
-# remove features with too low between-subject variance
 gm_maps_masked = nifti_masker.fit_transform(gray_matter_map_filenames)
-gm_maps_masked[:, gm_maps_masked.var(0) < 0.01] = 0.
-# final masking
-new_images = nifti_masker.inverse_transform(gm_maps_masked)
-gm_maps_masked = nifti_masker.fit_transform(new_images)
 n_samples, n_features = gm_maps_masked.shape
 print("%d samples, %d features" % (n_subjects, n_features))
 
@@ -80,7 +75,11 @@
 svr = SVR(kernel='linear')
 
 # Dimension reduction
-from sklearn.feature_selection import SelectKBest, f_regression
+from sklearn.feature_selection import VarianceThreshold, SelectKBest, \
+        f_regression
+
+# Remove features with too low between-subject variance
+variance_threshold = VarianceThreshold(threshold=.01)
 
 # Here we use a classical univariate feature selection based on F-test,
 # namely Anova.
@@ -90,7 +89,10 @@
 # we can plug them together in a *pipeline* that performs the two operations
 # successively:
 from sklearn.pipeline import Pipeline
-anova_svr = Pipeline([('anova', feature_selection), ('svr', svr)])
+anova_svr = Pipeline([
+            ('variance_threshold', variance_threshold),
+            ('anova', feature_selection),
+            ('svr', svr)])
 
 ### Fit and predict
 anova_svr.fit(gm_maps_masked, age)
@@ -101,6 +103,8 @@
 coef = svr.coef_
 # reverse feature selection
 coef = feature_selection.inverse_transform(coef)
+# reverse variance threshold
+coef = variance_threshold.inverse_transform(coef)
 # reverse masking
 weight_img = nifti_masker.inverse_transform(coef)
 
@@ -137,13 +141,14 @@
 
 # Statistical inference
 from nilearn.mass_univariate import permuted_ols
+data = variance_threshold.fit_transform(gm_maps_masked)
 neg_log_pvals, t_scores_original_data, _ = permuted_ols(
-    age, gm_maps_masked,  # + intercept as a covariate by default
-    n_perm=1000,  # 1,000 in the interest of time; 10000 would be better
+    age, data,  # + intercept as a covariate by default
+    n_perm=2000,  # 1,000 in the interest of time; 10000 would be better
     n_jobs=1)  # can be changed to use more CPUs
 signed_neg_log_pvals = neg_log_pvals * np.sign(t_scores_original_data)
 signed_neg_log_pvals_unmasked = nifti_masker.inverse_transform(
-    signed_neg_log_pvals)
+    variance_threshold.inverse_transform(signed_neg_log_pvals))
 
 # Show results
 threshold = -np.log10(0.1)  # 10% corrected

{"golden_diff": "diff --git a/examples/decoding/plot_oasis_vbm.py b/examples/decoding/plot_oasis_vbm.py\n--- a/examples/decoding/plot_oasis_vbm.py\n+++ b/examples/decoding/plot_oasis_vbm.py\n@@ -63,12 +63,7 @@\n     standardize=False,\n     smoothing_fwhm=2,\n     memory='nilearn_cache')  # cache options\n-# remove features with too low between-subject variance\n gm_maps_masked = nifti_masker.fit_transform(gray_matter_map_filenames)\n-gm_maps_masked[:, gm_maps_masked.var(0) < 0.01] = 0.\n-# final masking\n-new_images = nifti_masker.inverse_transform(gm_maps_masked)\n-gm_maps_masked = nifti_masker.fit_transform(new_images)\n n_samples, n_features = gm_maps_masked.shape\n print(\"%d samples, %d features\" % (n_subjects, n_features))\n \n@@ -80,7 +75,11 @@\n svr = SVR(kernel='linear')\n \n # Dimension reduction\n-from sklearn.feature_selection import SelectKBest, f_regression\n+from sklearn.feature_selection import VarianceThreshold, SelectKBest, \\\n+        f_regression\n+\n+# Remove features with too low between-subject variance\n+variance_threshold = VarianceThreshold(threshold=.01)\n \n # Here we use a classical univariate feature selection based on F-test,\n # namely Anova.\n@@ -90,7 +89,10 @@\n # we can plug them together in a *pipeline* that performs the two operations\n # successively:\n from sklearn.pipeline import Pipeline\n-anova_svr = Pipeline([('anova', feature_selection), ('svr', svr)])\n+anova_svr = Pipeline([\n+            ('variance_threshold', variance_threshold),\n+            ('anova', feature_selection),\n+            ('svr', svr)])\n \n ### Fit and predict\n anova_svr.fit(gm_maps_masked, age)\n@@ -101,6 +103,8 @@\n coef = svr.coef_\n # reverse feature selection\n coef = feature_selection.inverse_transform(coef)\n+# reverse variance threshold\n+coef = variance_threshold.inverse_transform(coef)\n # reverse masking\n weight_img = nifti_masker.inverse_transform(coef)\n \n@@ -137,13 +141,14 @@\n \n # Statistical inference\n from nilearn.mass_univariate import permuted_ols\n+data = variance_threshold.fit_transform(gm_maps_masked)\n neg_log_pvals, t_scores_original_data, _ = permuted_ols(\n-    age, gm_maps_masked,  # + intercept as a covariate by default\n-    n_perm=1000,  # 1,000 in the interest of time; 10000 would be better\n+    age, data,  # + intercept as a covariate by default\n+    n_perm=2000,  # 1,000 in the interest of time; 10000 would be better\n     n_jobs=1)  # can be changed to use more CPUs\n signed_neg_log_pvals = neg_log_pvals * np.sign(t_scores_original_data)\n signed_neg_log_pvals_unmasked = nifti_masker.inverse_transform(\n-    signed_neg_log_pvals)\n+    variance_threshold.inverse_transform(signed_neg_log_pvals))\n \n # Show results\n threshold = -np.log10(0.1)  # 10% corrected\n", "issue": "Misleading (wrong?) code in Oasis VBM example\nAs pointed out on http://neurostars.org/p/3388/, the Oasis VBM code applies the NiftiMasker twice, and thus ends up smoothing the data twice.\n\nIt does that to do screening of low variance features. I think that the right way to do this would be to use a sklearn.feature_selection.VarianceThreshold in the pipeline.\n\n", "before_files": [{"content": "\"\"\"\nVoxel-Based Morphometry on Oasis dataset\n========================================\n\nThis example uses Voxel-Based Morphometry (VBM) to study the relationship\nbetween aging and gray matter density.\n\nThe data come from the `OASIS <http://www.oasis-brains.org/>`_ project.\nIf you use it, you need to agree with the data usage agreement available\non the website.\n\nIt has been run through a standard VBM pipeline (using SPM8 and\nNewSegment) to create VBM maps, which we study here.\n\nPredictive modeling analysis: VBM bio-markers of aging?\n--------------------------------------------------------\n\nWe run a standard SVM-ANOVA nilearn pipeline to predict age from the VBM\ndata. We use only 100 subjects from the OASIS dataset to limit the memory\nusage.\n\nNote that for an actual predictive modeling study of aging, the study\nshould be ran on the full set of subjects. Also, parameters such as the\nsmoothing applied to the data and the number of features selected by the\nAnova step should be set by nested cross-validation, as they impact\nsignificantly the prediction score.\n\nBrain mapping with mass univariate\n-----------------------------------\n\nSVM weights are very noisy, partly because heavy smoothing is detrimental\nfor the prediction here. A standard analysis using mass-univariate GLM\n(here permuted to have exact correction for multiple comparisons) gives a\nmuch clearer view of the important regions.\n\n____\n\n\"\"\"\n# Authors: Elvis Dhomatob, <elvis.dohmatob@inria.fr>, Apr. 2014\n#          Virgile Fritsch, <virgile.fritsch@inria.fr>, Apr 2014\n#          Gael Varoquaux, Apr 2014\nimport numpy as np\nfrom scipy import linalg\nimport matplotlib.pyplot as plt\nfrom nilearn import datasets\nfrom nilearn.input_data import NiftiMasker\n\nn_subjects = 100  # more subjects requires more memory\n\n### Load Oasis dataset ########################################################\noasis_dataset = datasets.fetch_oasis_vbm(n_subjects=n_subjects)\ngray_matter_map_filenames = oasis_dataset.gray_matter_maps\nage = oasis_dataset.ext_vars['age'].astype(float)\n\n# print basic information on the dataset\nprint('First gray-matter anatomy image (3D) is located at: %s' %\n      oasis_dataset.gray_matter_maps[0])  # 3D data\nprint('First white-matter anatomy image (3D) is located at: %s' %\n      oasis_dataset.white_matter_maps[0])  # 3D data\n\n### Preprocess data ###########################################################\nnifti_masker = NiftiMasker(\n    standardize=False,\n    smoothing_fwhm=2,\n    memory='nilearn_cache')  # cache options\n# remove features with too low between-subject variance\ngm_maps_masked = nifti_masker.fit_transform(gray_matter_map_filenames)\ngm_maps_masked[:, gm_maps_masked.var(0) < 0.01] = 0.\n# final masking\nnew_images = nifti_masker.inverse_transform(gm_maps_masked)\ngm_maps_masked = nifti_masker.fit_transform(new_images)\nn_samples, n_features = gm_maps_masked.shape\nprint(\"%d samples, %d features\" % (n_subjects, n_features))\n\n### Prediction with SVR #######################################################\nprint(\"ANOVA + SVR\")\n# Define the prediction function to be used.\n# Here we use a Support Vector Classification, with a linear kernel\nfrom sklearn.svm import SVR\nsvr = SVR(kernel='linear')\n\n# Dimension reduction\nfrom sklearn.feature_selection import SelectKBest, f_regression\n\n# Here we use a classical univariate feature selection based on F-test,\n# namely Anova.\nfeature_selection = SelectKBest(f_regression, k=2000)\n\n# We have our predictor (SVR), our feature selection (SelectKBest), and now,\n# we can plug them together in a *pipeline* that performs the two operations\n# successively:\nfrom sklearn.pipeline import Pipeline\nanova_svr = Pipeline([('anova', feature_selection), ('svr', svr)])\n\n### Fit and predict\nanova_svr.fit(gm_maps_masked, age)\nage_pred = anova_svr.predict(gm_maps_masked)\n\n# Visualization\n# Look at the SVR's discriminating weights\ncoef = svr.coef_\n# reverse feature selection\ncoef = feature_selection.inverse_transform(coef)\n# reverse masking\nweight_img = nifti_masker.inverse_transform(coef)\n\n# Create the figure\nfrom nilearn.plotting import plot_stat_map, show\nbg_filename = gray_matter_map_filenames[0]\nz_slice = 0\nfrom nilearn.image.resampling import coord_transform\naffine = weight_img.get_affine()\n_, _, k_slice = coord_transform(0, 0, z_slice,\n                                linalg.inv(affine))\nk_slice = np.round(k_slice)\n\nfig = plt.figure(figsize=(5.5, 7.5), facecolor='k')\nweight_slice_data = weight_img.get_data()[..., k_slice, 0]\nvmax = max(-np.min(weight_slice_data), np.max(weight_slice_data)) * 0.5\ndisplay = plot_stat_map(weight_img, bg_img=bg_filename,\n                        display_mode='z', cut_coords=[z_slice],\n                        figure=fig, vmax=vmax)\ndisplay.title('SVM weights', y=1.2)\n\n# Measure accuracy with cross validation\nfrom sklearn.cross_validation import cross_val_score\ncv_scores = cross_val_score(anova_svr, gm_maps_masked, age)\n\n# Return the corresponding mean prediction accuracy\nprediction_accuracy = np.mean(cv_scores)\nprint(\"=== ANOVA ===\")\nprint(\"Prediction accuracy: %f\" % prediction_accuracy)\nprint(\"\")\n\n### Inference with massively univariate model #################################\nprint(\"Massively univariate model\")\n\n# Statistical inference\nfrom nilearn.mass_univariate import permuted_ols\nneg_log_pvals, t_scores_original_data, _ = permuted_ols(\n    age, gm_maps_masked,  # + intercept as a covariate by default\n    n_perm=1000,  # 1,000 in the interest of time; 10000 would be better\n    n_jobs=1)  # can be changed to use more CPUs\nsigned_neg_log_pvals = neg_log_pvals * np.sign(t_scores_original_data)\nsigned_neg_log_pvals_unmasked = nifti_masker.inverse_transform(\n    signed_neg_log_pvals)\n\n# Show results\nthreshold = -np.log10(0.1)  # 10% corrected\n\nfig = plt.figure(figsize=(5.5, 7.5), facecolor='k')\n\ndisplay = plot_stat_map(signed_neg_log_pvals_unmasked, bg_img=bg_filename,\n                        threshold=threshold, cmap=plt.cm.RdBu_r,\n                        display_mode='z', cut_coords=[z_slice],\n                        figure=fig)\ntitle = ('Negative $\\log_{10}$ p-values'\n         '\\n(Non-parametric + max-type correction)')\ndisplay.title(title, y=1.2)\n\nsigned_neg_log_pvals_slice_data = \\\n    signed_neg_log_pvals_unmasked.get_data()[..., k_slice, 0]\nn_detections = (np.abs(signed_neg_log_pvals_slice_data) > threshold).sum()\nprint('\\n%d detections' % n_detections)\n\nshow()\n", "path": "examples/decoding/plot_oasis_vbm.py"}], "after_files": [{"content": "\"\"\"\nVoxel-Based Morphometry on Oasis dataset\n========================================\n\nThis example uses Voxel-Based Morphometry (VBM) to study the relationship\nbetween aging and gray matter density.\n\nThe data come from the `OASIS <http://www.oasis-brains.org/>`_ project.\nIf you use it, you need to agree with the data usage agreement available\non the website.\n\nIt has been run through a standard VBM pipeline (using SPM8 and\nNewSegment) to create VBM maps, which we study here.\n\nPredictive modeling analysis: VBM bio-markers of aging?\n--------------------------------------------------------\n\nWe run a standard SVM-ANOVA nilearn pipeline to predict age from the VBM\ndata. We use only 100 subjects from the OASIS dataset to limit the memory\nusage.\n\nNote that for an actual predictive modeling study of aging, the study\nshould be ran on the full set of subjects. Also, parameters such as the\nsmoothing applied to the data and the number of features selected by the\nAnova step should be set by nested cross-validation, as they impact\nsignificantly the prediction score.\n\nBrain mapping with mass univariate\n-----------------------------------\n\nSVM weights are very noisy, partly because heavy smoothing is detrimental\nfor the prediction here. A standard analysis using mass-univariate GLM\n(here permuted to have exact correction for multiple comparisons) gives a\nmuch clearer view of the important regions.\n\n____\n\n\"\"\"\n# Authors: Elvis Dhomatob, <elvis.dohmatob@inria.fr>, Apr. 2014\n#          Virgile Fritsch, <virgile.fritsch@inria.fr>, Apr 2014\n#          Gael Varoquaux, Apr 2014\nimport numpy as np\nfrom scipy import linalg\nimport matplotlib.pyplot as plt\nfrom nilearn import datasets\nfrom nilearn.input_data import NiftiMasker\n\nn_subjects = 100  # more subjects requires more memory\n\n### Load Oasis dataset ########################################################\noasis_dataset = datasets.fetch_oasis_vbm(n_subjects=n_subjects)\ngray_matter_map_filenames = oasis_dataset.gray_matter_maps\nage = oasis_dataset.ext_vars['age'].astype(float)\n\n# print basic information on the dataset\nprint('First gray-matter anatomy image (3D) is located at: %s' %\n      oasis_dataset.gray_matter_maps[0])  # 3D data\nprint('First white-matter anatomy image (3D) is located at: %s' %\n      oasis_dataset.white_matter_maps[0])  # 3D data\n\n### Preprocess data ###########################################################\nnifti_masker = NiftiMasker(\n    standardize=False,\n    smoothing_fwhm=2,\n    memory='nilearn_cache')  # cache options\ngm_maps_masked = nifti_masker.fit_transform(gray_matter_map_filenames)\nn_samples, n_features = gm_maps_masked.shape\nprint(\"%d samples, %d features\" % (n_subjects, n_features))\n\n### Prediction with SVR #######################################################\nprint(\"ANOVA + SVR\")\n# Define the prediction function to be used.\n# Here we use a Support Vector Classification, with a linear kernel\nfrom sklearn.svm import SVR\nsvr = SVR(kernel='linear')\n\n# Dimension reduction\nfrom sklearn.feature_selection import VarianceThreshold, SelectKBest, \\\n        f_regression\n\n# Remove features with too low between-subject variance\nvariance_threshold = VarianceThreshold(threshold=.01)\n\n# Here we use a classical univariate feature selection based on F-test,\n# namely Anova.\nfeature_selection = SelectKBest(f_regression, k=2000)\n\n# We have our predictor (SVR), our feature selection (SelectKBest), and now,\n# we can plug them together in a *pipeline* that performs the two operations\n# successively:\nfrom sklearn.pipeline import Pipeline\nanova_svr = Pipeline([\n            ('variance_threshold', variance_threshold),\n            ('anova', feature_selection),\n            ('svr', svr)])\n\n### Fit and predict\nanova_svr.fit(gm_maps_masked, age)\nage_pred = anova_svr.predict(gm_maps_masked)\n\n# Visualization\n# Look at the SVR's discriminating weights\ncoef = svr.coef_\n# reverse feature selection\ncoef = feature_selection.inverse_transform(coef)\n# reverse variance threshold\ncoef = variance_threshold.inverse_transform(coef)\n# reverse masking\nweight_img = nifti_masker.inverse_transform(coef)\n\n# Create the figure\nfrom nilearn.plotting import plot_stat_map, show\nbg_filename = gray_matter_map_filenames[0]\nz_slice = 0\nfrom nilearn.image.resampling import coord_transform\naffine = weight_img.get_affine()\n_, _, k_slice = coord_transform(0, 0, z_slice,\n                                linalg.inv(affine))\nk_slice = np.round(k_slice)\n\nfig = plt.figure(figsize=(5.5, 7.5), facecolor='k')\nweight_slice_data = weight_img.get_data()[..., k_slice, 0]\nvmax = max(-np.min(weight_slice_data), np.max(weight_slice_data)) * 0.5\ndisplay = plot_stat_map(weight_img, bg_img=bg_filename,\n                        display_mode='z', cut_coords=[z_slice],\n                        figure=fig, vmax=vmax)\ndisplay.title('SVM weights', y=1.2)\n\n# Measure accuracy with cross validation\nfrom sklearn.cross_validation import cross_val_score\ncv_scores = cross_val_score(anova_svr, gm_maps_masked, age)\n\n# Return the corresponding mean prediction accuracy\nprediction_accuracy = np.mean(cv_scores)\nprint(\"=== ANOVA ===\")\nprint(\"Prediction accuracy: %f\" % prediction_accuracy)\nprint(\"\")\n\n### Inference with massively univariate model #################################\nprint(\"Massively univariate model\")\n\n# Statistical inference\nfrom nilearn.mass_univariate import permuted_ols\ndata = variance_threshold.fit_transform(gm_maps_masked)\nneg_log_pvals, t_scores_original_data, _ = permuted_ols(\n    age, data,  # + intercept as a covariate by default\n    n_perm=2000,  # 1,000 in the interest of time; 10000 would be better\n    n_jobs=1)  # can be changed to use more CPUs\nsigned_neg_log_pvals = neg_log_pvals * np.sign(t_scores_original_data)\nsigned_neg_log_pvals_unmasked = nifti_masker.inverse_transform(\n    variance_threshold.inverse_transform(signed_neg_log_pvals))\n\n# Show results\nthreshold = -np.log10(0.1)  # 10% corrected\n\nfig = plt.figure(figsize=(5.5, 7.5), facecolor='k')\n\ndisplay = plot_stat_map(signed_neg_log_pvals_unmasked, bg_img=bg_filename,\n                        threshold=threshold, cmap=plt.cm.RdBu_r,\n                        display_mode='z', cut_coords=[z_slice],\n                        figure=fig)\ntitle = ('Negative $\\log_{10}$ p-values'\n         '\\n(Non-parametric + max-type correction)')\ndisplay.title(title, y=1.2)\n\nsigned_neg_log_pvals_slice_data = \\\n    signed_neg_log_pvals_unmasked.get_data()[..., k_slice, 0]\nn_detections = (np.abs(signed_neg_log_pvals_slice_data) > threshold).sum()\nprint('\\n%d detections' % n_detections)\n\nshow()\n", "path": "examples/decoding/plot_oasis_vbm.py"}]}