Ich habe versucht, fMRI-Daten mit Python zu analysieren (Einführung in die Dekodierung von Gehirninformationen)
Qiita hat kein Analysebeispiel für fMRI-Daten ... !!
Neulich dachte ich: "Es gibt verschiedene Artikel über Qiita ... Ich denke, es gibt nicht viele Beispiele für die Analyse von ** fMRI-Daten ** mit Python oder maschinellem Lernen." Daher suchte ich in Qiita nach Beispielen für die tatsächliche Analyse von fMRI-Daten von fMRI.
Das? Nicht wirklich ...? !! Es gibt Kommentarartikel in den Bereichen Werkzeuge und Neurowissenschaften, aber ich konnte keine Analyse unter Verwendung der tatsächlichen fMRI-Daten finden.
Also ** dieses Mal werden wir maschinelles Lernen unter Verwendung von fMRI-Daten ** analysieren. Obwohl es ein Gehirnbild ist, ist es doch ein dreidimensionales Bild. Eine Querschnittsansicht eines zweidimensionalen Bildes, das durch Halbieren des Gehirns erhalten wird, wird durch Überlappen gebildet, und jedes Pixel enthält Informationen in einer Einheit, die als "** boxel **" bezeichnet wird. Ich bin nur da. Wenn es um 4D-Gehirnbilder geht, enthält jede dreidimensionale Boxel-Information nur Zeitreiheninformationen (beim Scannen mit fMRI).
Gehirnbilder werden häufig in einem speziellen Format veröffentlicht, das als NIfTY-Format bezeichnet wird und gelesen werden muss. Dieses Mal verwenden wir anstelle von Matlab, das angeblich häufig von Studenten und Forschern verwendet wird, die ** Nilearn ** -Bibliothek von Python, um Gehirnbilder zu lesen und zu analysieren (denken Sie daran, dass jeder sie analysieren kann). Aus diesem Grund habe ich aufgehört, Matlab zu bezahlen. Die Analyse des maschinellen Lernens selbst wird von der SVM ** scikit-learn ** gelernt und vorhergesagt, die jeder, der maschinelles Lernen durchgeführt hat, verwendet oder davon gehört hat.
Obwohl es sich um einen Ansatz handelt, der maschinelles Lernen verwendet, werden wir dieses Mal fMRI-Daten verwenden, die eine bestimmte Gehirnregion maskieren, um eine Analyse durchzuführen, die vorhersagt, "was fMRI zu diesem Zeitpunkt aufgrund der Gehirnaktivität gesehen hat" * *. In der Gehirnanalysebranche scheint es als ** Dekodierung von Gehirninformationen ** bezeichnet zu werden, aber es ist nur maschinelles Lernen (Entschuldigung für die Leute in der Branche). In der Vergangenheit verwendeten wir einen Ansatz wie die Codierung, bei dem das in einer experimentellen Aufgabe durchgeführte Verhalten als erklärende Variable verwendet wurde und die Gehirnregion, die statistisch signifikant reagierte, wenn sie von fMRI durchgeführt wurde, als unabhängige Variable erkannt wurde. Ich denke, es kam von der Tatsache, dass es Dekodierung genannt wurde, weil es den umgekehrten Ansatz verfolgte (das heißt, anstatt die Gehirnregion zu betrachten, die von der Aufgabe statistisch signifikant angeregt wurde, aber die Aktivität einer Gehirnregion. Vorhersage der Aufgabe von).
Über die tatsächlich verwendeten Daten
Dieses Mal analysierte ich es anhand des Tutorials zur Decodierungsanalyse von Nilearn: Maschinelles Lernen für Neuro-Imaging in Python, das zum Lesen von fMRI-Daten usw. verwendet wird. ** Das Haxby 2001-Experiment ** Verwenden Sie das in veröffentlichte Gehirnbild.
A introduction tutorial to fMRI decoding
Wenn Sie dem obigen Tutorial folgen, heißt es, dass Sie mit dem folgenden Code problemlos einen Datensatz von Probanden aus dem Haxby 2001-Experiment erfassen können.
from nilearn import datasets
#Standardmäßig werden Daten für das zweite Thema abgerufen
haxby_dataset = datasets.fetch_haxby()
Es scheint jedoch, dass es ein Problem mit der NIRTC-Datenbank gibt, die die Daten bereitstellt, aber ich konnte sie nicht herunterladen. Daher habe ich sie direkt heruntergeladen und lokal gespeichert. Sie können subj2-2010.01.14.tar.gz, das die Daten des zweiten Betreffs enthält, von dieser URL herunterladen. Wenn Sie die Daten analysieren möchten, laden Sie sie bitte von hier herunter (auch Daten anderer Betreffnummern). macht nichts).
http://data.pymvpa.org/datasets/haxby2001/
In diesem Experiment werden den Probanden verschiedene Kategorien visueller Reize (wie Scheren, menschliche Gesichter und Katzen) in einem fMRT-Scanner präsentiert. In diesem Tutorial sagt die im Bereich des ventralen kortikalen visuellen Trakts (der an der Darstellung von Farbe und Form beteiligt ist) aufgezeichnete fMRI-Gehirnaktivität voraus, welche Kategorie das Subjekt mit fMRI betrachtet.
Analyseverfahren
- Im Quellcode wurde der Name des heruntergeladenen Ordners in subj2 geändert.
Laden und visualisieren Sie Bilder
Lesen wir zunächst die 4D-fMRI-Daten. Wie oben erwähnt, ist dies ein Format, das Zeitreiheninformationen enthält, wenn eine experimentelle Aufgabe in einem Gehirnbild ausgeführt wird, das Boxeldaten (dreidimensionale Daten) sind.
Da es sich jedoch um 4D-fMRI-Daten handelt, können sie nicht durch Verarbeitung entsprechender 3D-fMRI-Daten (ohne Zeitreiheninformationen) wie nilearn.plotting.plot_epi visualisiert werden. Daher werden die visualisierten 4D-fMRI-Daten gemittelt und in eine 3D-fMRI-Daten umgewandelt. Dies kann durch Importieren von "from nilearn.image import mean_img" erreicht werden.
fmri_filename = "~/Desktop/nilearn/subj2/bold.nii.gz"
from nilearn import plotting
from nilearn.image import mean_img
plotting.view_img(mean_img(fmri_filename), threshold=None)
Ergebnis
Zunächst konnte ich das Gehirnbild visualisieren. Es ist großartig, in wenigen Zeilen visualisieren zu können.
Extrahieren Sie Funktionen von fMRI-Daten und konvertieren Sie sie in Daten im Numpy-Format
Vielleicht waren einige Leute erleichtert, das Wort numpy zu sehen. Hier verwenden wir eine Funktion namens "nilearn.input_data.NiftiMasker", um die fMRI-Daten des ventralen kortikalen visuellen Trakts zu maskieren und in das Numpy-Format zu konvertieren.
Lassen Sie uns zunächst die fMRT-Daten nur des ventralen kortikalen visuellen Trakts visualisieren.
mask_filename = "~/Desktop/nilearn/subj2/mask4_vt.nii.gz"
anat_filename = "~/Desktop/nilearn/subj2/anat.nii.gz"
#Versuchen Sie, den Bereich der Maske vor dem Hintergrund des anatomischen Bilds des Motivs anzuzeigen
plotting.plot_roi(mask_filename, bg_img=anat_filename,
cmap='Paired')
Ergebnis
Als nächstes erstellen wir eine Maske. Ich möchte die Daten standardisieren, um die Genauigkeit des maschinellen Lernens zu verbessern. Die Maske wird in Form eines 2D-Arrays extrahiert, wenn sie für das maschinelle Lernen mit nilearn bereit ist.
from nilearn.input_data import NiftiMasker
masker = NiftiMasker(mask_img=mask_filename, standardize=True) #Standardisieren Sie die Maske
fmri_masked = masker.fit_transform(fmri_filename) #Legen Sie das BOLD-Signal auf die standardisierte Maske an
Übrigens können Sie durch Ausführen von masker.generate_report () das NIfTY-Bild der Maske anzeigen und verschiedene Parameter usw. überprüfen.
Sie haben jetzt die maskierte Datenmatrix in der Form fmri_masked als Numpy-Format extrahiert. Ich werde es tatsächlich versuchen.
# fmri_maskiert ist in Form eines numpy Arrays.
print(fmri_masked)
#Die Größe dieser Daten ist die Anzahl der Boxen und Zeitreihenpunkte(Anzahl der Scans)Es besteht aus der Gesamtzahl von.
# shape[0]Ist die Anzahl der Scans, shape[1]Ist die Anzahl der Boxen
print(fmri_masked.shape)
Ausgabeergebnis
[[ 7.6757896e-01 2.3108697e+00 -2.0519443e-01 ... -1.0261143e+00
8.7993503e-02 2.0720518e+00]
[ 5.5640817e-01 1.6833434e+00 -2.4644937e-01 ... -7.0238107e-01
-3.4570047e-01 2.0341001e+00]
[ 7.6757896e-01 1.9186659e+00 1.0802225e-03 ... -9.9374104e-01
-2.7630943e-01 2.1479552e+00]
...
[-4.2905563e-01 -1.6896105e+00 -7.4150854e-01 ... -1.5440876e+00
1.8054217e+00 -1.6709718e-01]
[-1.4749455e-01 -1.8072717e+00 -2.4644937e-01 ... -1.7707009e+00
1.5452052e+00 7.8169477e-01]
[-2.1788482e-01 -1.4542881e+00 1.0802225e-03 ... -1.6412076e+00
1.2676411e+00 8.9554977e-01]]
(1452, 464)
Wie durch "fmri_masked.shape" bestätigt, speichert 464 die Anzahl der Boxzellen in der maskierten Gehirnregion und 1452 speichert Zeitreiheninformationen (fMRI wird mit Zeitauflösung in Einheiten von TR erhalten und fMRI-Daten werden erfasst. Normalerweise wird es in Intervallen von etwa 1 bis 2,5 Sekunden eingestellt, aber die Dekodierung analysiert anhand von Informationen zu einer bestimmten Gehirnregion und Neurofeedback, die die Gehirnaktivität in Echtzeit messen und das Subjekt trainieren Bei der aufgerufenen Methode denke ich, dass oft eine ziemlich kurze TR eingestellt wurde.
Wenn Sie als Test die Zeitreiheninformationen der ersten drei Boxen herausnehmen, sieht dies wie folgt aus.
#Zeigen Sie die ersten drei aus der Matrix ausgewählten Boxen an.
import matplotlib.pyplot as plt
plt.plot(fmri_masked[5:150, :3])
plt.title('Voxel Time Series')
plt.xlabel('Scan number')
plt.ylabel('Normalized signal')
plt.tight_layout()
Übrigens mögen sich scharfe Leute fragen, warum der Teil "fmri_masked [5: 150,: 3]" nicht die Form von "fmri_masked [: 150,: 3]" hat. Dies wird als ** Dummy-Scan ** bezeichnet, da der erste Scan, der mit der Abbildung von fMRI beginnt, einen hohen Signalwert aufweist und in der tatsächlichen Analyse häufig nicht verwendet wird. Der Signalwert ist hoch, da der T1-Effekt nicht stabil ist (Details finden Sie auf anderen Websites bei Google).
Erfassen Sie Aktionsbeschriftungen während experimenteller Aufgaben
Nachdem wir die fMRI-Daten des ventralen kortikalen visuellen Trakts im Numpy-Format extrahiert haben, werden wir diese Daten für maschinelles Lernen verwenden. ** Dieses Mal lernen wir mit einem Lehrer, daher benötigen wir zusätzlich zu den Daten ein Etikett **. Daher werden wir daran arbeiten, die Informationen der Kategorie, die das Subjekt zum Zeitpunkt des Experiments betrachtete, in den heruntergeladenen Daten als Etikett zu extrahieren.
import pandas as pd
#Lesen Sie die Informationen zum Aktionsergebnis. Es gibt so viele Etiketten wie Scans
behavioral = pd.read_csv('~/Desktop/nilearn/subj2/labels.txt', delimiter=' ')
Ausgabeergebnis
labels chunks
0 rest 0
1 rest 0
2 rest 0
3 rest 0
4 rest 0
5 rest 0
6 scissors 0
7 scissors 0
8 scissors 0
9 scissors 0
10 scissors 0
11 scissors 0
12 scissors 0
13 scissors 0
14 scissors 0
... ... ...
1426 cat 11
1427 cat 11
1428 cat 11
1429 cat 11
1430 cat 11
1431 cat 11
1432 rest 11
1433 rest 11
1434 rest 11
1435 rest 11
1436 rest 11
1437 scissors 11
1438 scissors 11
1439 scissors 11
1440 scissors 11
1441 scissors 11
1442 scissors 11
1443 scissors 11
1444 scissors 11
1445 scissors 11
1446 rest 11
1447 rest 11
1448 rest 11
1449 rest 11
1450 rest 11
1451 rest 11
1452 rows × 2 columns
Da es sich bei dieser Aufgabe um eine visuelle Erkennungsaufgabe handelt, wie oben erwähnt, zeigt das Etikett die experimentellen Bedingungen (Kategorie des Objekts, das dem Subjekt gezeigt wird). Chunks stellen die Anzahl der Sitzungen dar (denn wenn Sie das Subjekt die Aufgabe für immer erledigen lassen, ist die Belastung des Körpers nicht ungerade. Normalerweise erledigen Sie ähnliche Aufgaben mehrmals, während Sie Ihre körperliche Verfassung in den Sitzungspausen überprüfen. Ich werde fortsetzen).
conditions = behavioral['labels']
Ausgabeergebnis
0 rest
1 rest
2 rest
3 rest
4 rest
5 rest
6 scissors
7 scissors
8 scissors
9 scissors
10 scissors
11 scissors
12 scissors
13 scissors
14 scissors
15 rest
16 rest
17 rest
18 rest
19 rest
20 rest
21 face
22 face
...
1431 cat
1432 rest
1433 rest
1434 rest
1435 rest
1436 rest
1437 scissors
1438 scissors
1439 scissors
1440 scissors
1441 scissors
1442 scissors
1443 scissors
1444 scissors
1445 scissors
1446 rest
1447 rest
1448 rest
1449 rest
1450 rest
1451 rest
Name: labels, Length: 1452, dtype: object
Dieses Mal werden wir nach dem Tutorial nur die experimentellen Bedingungen für Katzen und Gesichter verwenden. Es gibt viele andere Labels. Wenn Sie also interessiert sind, versuchen Sie es bitte in anderen Kategorien.
#Holen Sie sich das Boxel, das der Anzahl der Scans unter den experimentellen Bedingungen entspricht, bei denen Gesicht und Katze dem maskierten Gehirnbild präsentiert werden.
condition_mask = conditions.isin(['face','cat'])
fmri_masked = fmri_masked[condition_mask]
print(fmri_masked.shape)
#Wenden Sie die gleichen Bedingungen auf Etiketten an
conditions = conditions[condition_mask]
print(conditions.shape)
Dekodierung mit SVM
Dieses Mal verwenden wir die Toolbox für maschinelles Lernen von scikit-learn für die abgerufenen fmri_masked-Daten. Ich benutze den linearen SVM-Kernel als Decoder. Da es sich um Scikit-Lernen handelt, können Sie alles vom Lernen bis zur Vorhersage in nur wenigen Zeilen erledigen.
from sklearn.svm import SVC
svc = SVC(kernel='linear')
#Lernen und Vorhersagen treffen
svc.fit(fmri_masked, conditions)
prediction = svc.predict(fmri_masked)
Dies ist jedoch für die Analyse völlig nutzlos, sodass wir (offensichtlich) eine Kreuzvalidierung durchführen werden.
Führen Sie eine Kreuzvalidierung durch, um die vorhergesagte Punktzahl zu berechnen
Führen Sie eine KFold-Kreuzungsüberprüfung durch
from sklearn.model_selection import KFold
import numpy as np
cv = KFold(n_splits=5)
scores = []
for train, test in cv.split(X=fmri_masked):
svc.fit(fmri_masked[train], conditions.values[train])
prediction = svc.predict(fmri_masked[test])
scores.append((prediction == conditions.values[test]).sum()
/ float(len(conditions.values[test])))
#Durchschnittliche Vorhersageleistung
print(np.mean(scores))
Ausgabeergebnis:
0.7628964059196617
Es ist auch möglich, die Funktion cross_val_score zu verwenden, um eine verteilte Verarbeitung unter Verwendung mehrerer CPUs auf der Maschine durchzuführen, um die Auswertungsverarbeitung für jede Teilmenge zu berechnen. (Sie kann in Form von n_jobs = n angegeben werden. Wenn -1 festgelegt ist, kann die Parallelverarbeitung mit allen auf dem Computer verfügbaren CPUs durchgeführt werden.)
from sklearn.model_selection import cross_val_score
cv_score = cross_val_score(svc, fmri_masked, conditions)
print(cv_score)
Ausgabeergebnis:
[0.86363636 0.76744186 0.74418605 0.69767442 0.81395349]
** Dies reicht jedoch noch nicht zur Überprüfung aus. Der Grund dafür ist, dass wir die Leistung nicht unter Berücksichtigung der Auswirkungen jeder Sitzung bewertet haben. ** ** **
Bewerten Sie die Vorhersagegenauigkeit für jede Sitzung
Die fMRI-Daten werden sitzungsweise erfasst, und das Rauschen ist für jede bestimmte Sitzung selbstkorreliert. Wenn Sie eine Kreuzvalidierung unter Berücksichtigung dieser Effekte durchführen, müssen Sie daher sitzungsübergreifende Vorhersagen treffen.
from sklearn.model_selection import LeaveOneGroupOut
cv = LeaveOneGroupOut()
cv_score = cross_val_score(svc,
fmri_masked,
conditions,
cv=cv,
groups=session_label,
)
#Erhalten Sie die Vorhersagegenauigkeit in den Verifizierungsdaten für jede Sitzung
print(cv_score)
Ausgabeergebnis:
[0.55555556 1. 0.66666667 0.66666667 0.77777778 0.72222222
0.88888889 0.38888889 0.66666667 0.5 0.77777778 0.66666667]
Schätzen Sie das Gewicht des trainierten Modells
Schließlich wird das Gewicht des Modells geschätzt und angezeigt. Dazu konvertieren wir zuerst die Gewichte in NIfTY-Bilder.
coef_ = svc.coef_
#Es liegt in Form eines Numpy-Arrays mit einem Koeffizienten pro Boxel vor.(Gewicht)haben
print(coef_.shape)
Ausgabeergebnis:
(1, 464)
coef_img = masker.inverse_transform(coef_)
print(coef_img)
Ausgabeergebnis:
<class 'nibabel.nifti1.Nifti1Image'>
data shape (40, 64, 64, 1)
affine:
[[ -3.5 0. 0. 68.25 ]
[ 0. 3.75 0. -118.125]
[ 0. 0. 3.75 -118.125]
[ 0. 0. 0. 1. ]]
metadata:
<class 'nibabel.nifti1.Nifti1Header'> object, endian='<'
sizeof_hdr : 348
data_type : b''
db_name : b''
extents : 0
session_error : 0
regular : b''
dim_info : 0
dim : [ 4 40 64 64 1 1 1 1]
intent_p1 : 0.0
intent_p2 : 0.0
intent_p3 : 0.0
intent_code : none
datatype : float64
bitpix : 64
slice_start : 0
pixdim : [-1. 3.5 3.75 3.75 1. 1. 1. 1. ]
vox_offset : 0.0
scl_slope : nan
scl_inter : nan
slice_end : 0
slice_code : unknown
xyzt_units : 0
cal_max : 0.0
cal_min : 0.0
slice_duration : 0.0
toffset : 0.0
glmax : 0
glmin : 0
descrip : b''
aux_file : b''
qform_code : unknown
sform_code : aligned
quatern_b : 0.0
quatern_c : 1.0
quatern_d : 0.0
qoffset_x : 68.25
qoffset_y : -118.125
qoffset_z : -118.125
srow_x : [-3.5 0. 0. 68.25]
srow_y : [ 0. 3.75 0. -118.125]
srow_z : [ 0. 0. 3.75 -118.125]
intent_name : b''
magic : b'n+1'
Speichern Sie das erstellte Bild im NIfTY-Format im Format nii.gz.
# nii.Speichern Sie im gz-Dateiformat
coef_img.to_filename('haxby_svc_weights.nii.gz')
Lassen Sie uns abschließend den Gewichtskoeffizienten als Gehirnbild visualisieren.
#Zeichnen Sie tatsächlich das Gewicht von SVM
from nilearn.plotting import plot_stat_map, show
plot_stat_map(coef_img, bg_img=anat_filename,
title="SVM weights", display_mode='yx')
show()
Mit dem oben Gesagten haben wir eine Reihe von fMRI-Datenanalysen unter Verwendung von maschinellem Lernen durchgeführt.
Zusammenfassung
Ich hatte vor, diese Analyse in meinem persönlichen Blog zu veröffentlichen, konnte aber niemanden finden, der fMRI-Daten mit Qiita analysierte. Deshalb habe ich sie hier veröffentlicht (einzelne Blogs werden in meinem Profil veröffentlicht).
Bild von fMRI-Daten ... Es mag schwierig erscheinen, es nur zu hören, aber am Ende ist es nur mehr Rauschen als ein allgemeines Bild, eine weitere Dimension und das Entfernen des Teils, der spezielle Formate verarbeitet Bildanalyse. Eines der Ziele ist es, in der Lage zu sein, verrückte Daten von jedermann zu analysieren. Ich veröffentliche sie über meinen persönlichen Blog und Qiita. Ich hoffe, dass Sie daran interessiert sind, Daten zu analysieren, die Sie durch diese Aktivität selten berühren.
Wenn es Ihnen nichts ausmacht, LG ... Vielen Dank für Ihre Mitarbeit! Lol
Eine Notiz namens Nachtrag
Dieses Mal werden wir anhand von rohen fMRI-Zeitreihendaten analysieren. Es ist jedoch erforderlich, die ursprünglich als HRF-Funktion bezeichnete Funktion zu falten und von der Aufgabenpräsentation (Nervenaktivität) in die Blutflussdynamik (BOLD-Reaktion) umzuwandeln. Um eine Test-für-Test-Aktivitätskarte (Beta-Karte) zu erhalten, ist es außerdem erforderlich, unter Verwendung eines allgemeinen linearen Modells auf der ersten Ebene zu passen. Dies ermöglicht es, Vorhersagen zwischen verwandten Bedingungen für jeden Versuch zu treffen. Tatsächlich sind verschiedene Verarbeitungsprozesse für die eigentliche fMRI-Datenanalyse erforderlich, aber diesmal handelt es sich um eine Einführung, die darauf abzielt, maschinelles Lernen auf fMRI-Daten anzuwenden. Daher habe ich mich entschlossen, die oben erwähnte Verarbeitung (tatsächlich) durchzuführen. Wenn Sie so viele Schritte unternehmen, passt es nicht in einen Artikel (Schweiß).
Bitte beachten Sie (ich möchte etwas ernsthaft analysieren, wenn es eine Gelegenheit gibt)
Recommended Posts