DE4326044C2

DE4326044C2 - Verfahren zur zeitgemittelten Lokalisierung von mit starkem Rauschen überlagerten elektrophysiologischen Aktivitäten

Info

Publication number: DE4326044C2
Application number: DE4326044A
Authority: DE
Inventors: Bernhard Dr Scholz; Sekihara Dr Kensuke
Original assignee: Siemens AG; Hitachi Ltd
Current assignee: Siemens AG; Hitachi Ltd
Priority date: 1992-10-19
Filing date: 1993-08-03
Publication date: 1995-08-17
Anticipated expiration: 2013-08-04
Also published as: DE4326044A1

Description

Die Erfindung betrifft ein Verfahren zur zeitgemittelten Lokalisierung von elektrophysiologischen Aktivitäten, die in einem Lebewesen innerhalb eines Untersuchungsgebiets auftre ten und die ein Magnetfeld erzeugen, das an Meßpunkten außer halb des Untersuchungsgebiets mit einem Vielkanalmeßsystem gemessen wird.

Ein Vielkanalmeßsystem, mit dem das eingangs genannte Magnet feld gemessen werden kann, ist aus der EP-A-0 359 864 bekannt. Das Vielkanalmeßsystem wird auch als biomagnetisches Meßsystem bezeichnet, mit dem sehr schwache Magnetfelder, die von im Innern eines Lebewesens ablaufenden elektrophysiolo gischen Aktivitäten erzeugt werden, gemessen werden können. Das Vielkanalmeßsystem umfaßt eine Vielkanal-Gradiometer anordnung, die mit einer Vielkanal-SQUID-Anordnung gekoppelt ist.

Mit dem biomagnetischen Meßsystem können sowohl Magnetoenze phalogramme (MEG) wie auch Magnetokardiogramme (MKG) gemessen werden. Das Hauptziel für die Auswertung der MEG- oder MKG- Aufzeichnungen ist eine dreidimensionale nicht-invasive Loka lisierung von Quellen elektrophysiologischer Aktivitäten.

In dem Artikel von Jukka Sarvas: "Basic mathematical and electromagnetic concepts of the biomagnetic inverse problem", erschienen in Phys. Med. Biol., 1987, Vol. 32, No. 1, pp. 11-22, ist in Kapitel 7 eine lineare Lösung für das magnetische inverse Problem, aus gemessenen Feldwerten den Ort und die Größe von die Felder erzeugenden Stromdichten zu bestimmen, angegeben. Die Stromdichten stellen ein Äquivalent der elek trophysiologischen Aktivität dar. Dort ist auch auf weiter gehende Ausführungen von R. J. Ilmoniemi und M. S. Hämäläinen in Helsinki University of Technology Report TKK-F-A559, 1984, hingewiesen. Danach wird ein Modell der elektrophysio logischen Aktivitäten das Untersuchungsgebiet in Zellen auf geteilt, in denen jeweils eine Stromdichte angenommen wird, die sich aus maximal drei Komponenten zusammensetzt. Dabei ist die Gesamtzahl der Komponenten größer als die Anzahl der Meßpunkte. Zwischen den Stromdichten und den Meßwerten an dem Meßorten besteht eine lineare Beziehung, die durch eine Lead- Field-Matrix beschrieben ist, die nur von der relativen Lage der Zellen zu den Meßorten abhängig ist. Die Lead-Field-Matrix stellt somit für jeden am Meßort sich befindenden Meß wertaufnehmer des Vielkanalmeßsystems ein Empfindlichkeits muster über das gesamte Untersuchungsgebiet dar. Das Unter suchungsgebiet wird auch als Rekonstruktionsgebiet bezeich net. Zur Rekonstruktion der bioelektrischen Stromdichtever teilung wird die Lead-Field-Matrix einer Pseudo-Inversion nach Moore-Penrose unterworfen. Diese invertierte Lead-Field- Matrix wird mit den an den Meßpunkten aufgenommenen Meßwerten multipliziert, um einen Minimum-Norm-Schätzwert der Stromdichte verteilung im Untersuchungsgebiet zu erhalten.

Biomagnetische Messungen haben jedoch oft ein niedriges Signal-Rausch-Verhältnis. Das Rauschen wird zum einen von externen Rauschquellen in der Umgebung und zum anderen von den SQUID-Gradiometern und der damit verbundenen Elektronik ver ursacht. Wenn das Signal-Rausch-Verhältnis der gemessenen Daten sehr niedrig ist, kann die Rekonstruktion über die nach Moore-Penrose gebildete Pseudo-Inversion der Lead-Field-Matrix keine aussagekräfigen Resultate liefern. Eine mögliche Lösung bei stark verrauschten Meßwerten besteht darin, die Meßwerte von sich wiederholenden Ereignissen zu mitteln. Bei einzelnen spontanen elektrophysiologischen Aktivitäten ist eine derartige Mittelung jedoch nicht möglich. Bei spontanen Aktivitäten müssen andere Verfahren eingesetzt werden, um Störungen im Meßsignal zu kompensieren.

Ein räumlich konstantes Störsignal wird nach einem aus der DE-OS 41 18 126 bekannten Verfahren eliminiert, indem von den Meßsignalen jedes Meßorts ein gewichtetes Summensignal sub trahiert wird, das aus den Meßsignalen zumindest einer Gruppe von Meßorten gebildet wird.

Aus der EP-A-0 477 434 ist bekannt, ein interessierendes Si gnal, das verrauscht ist und zusätzlich von anderen Signalen überlagert ist, zu erkennen und zu isolieren. Dazu werden die an den Meßorten gemessenen Feldwerte gewichtet, um einen vir tuellen Sensor zu definieren, der zur Messung eines bestimm ten Ortes der Aktivität optimiert ist. Die Gewichtskoeffi zienten werden aus der Kovarianzmatrix der Meßwerte an den Meßorten oder gemäß einem mathematischen Modell bestimmt, das Informationen über die Aktivität, das Untersuchungsgebiet und den Orten und der Orientierung der Aktivitäten und der Sensoren des Vielkanalmeßsystems benutzt.

Weitere Möglichkeiten das Signal-Rausch-Verhältnis zu verbes sern ist aus Samuel J. Williamson et al.: Advances in Bioma gnetism, Plenum Press, New York ans London, 1989, Proc. of the Seventh International Conference on Biomagnetism, August 13-18, 1989, New York, pp. 721-736, bekannt. Danach besteht eine Möglichkeit darin, mit einem Referenzsensor das Störsi gnal zu erfassen und ggf. gewichtet von den Meßwerten abzu ziehen. Untersuchungen über verschiedene Verfahren zur Be stimmung der Gewichtsfaktoren sind ebenfalls beschrieben. Thermisches magnetisches Rauschen von elektrisch leitfähigen Materialien bilden auch eine Störquelle. So sind dort Studien beschrieben, worin das thermische magnetische Rauschen ein schließlich der räumlichen Korrelation berechnet wird, das von einem Dewar erzeugt wird.

Bei im wesentlichen ortsfesten elektrophysiologischen Akti vitäten genügt es oft, nur eine über ein Meßintervall zeit lich gemittelte Rekonstruktion durchzuführen.

Der Erfindung liegt nun die Aufgabe zugrunde, ein Verfahren anzugeben, mit dem der Ort von zeitgemittelten elektrophy siologischen Aktivitäten innerhalb eines Untersuchungsge biets bestimmt werden kann, wenn die außerhalb des Unter suchungsgebiets an Meßorten ermittelten Meßwerte des von den Aktivitäten erzeugten Magnetfeldes von einem starken Rauschen überlagert sind.

Die Aufgabe wird gelöst durch ein Verfahren mit den Schritten:

a) in einem Modell der elektrophysiologischen Aktivitäten im Untersuchungsgebiet wird das Untersuchungsgebiet in Zellen aufgeteilt, in denen jeweils eine Stromdichte vorhanden ist, die sich aus maximal drei Komponenten zusammensetzt, wobei die Gesamtzahl der Komponenten größer ist als die Anzahl der Meßpunkte und wobei zwi schen den Stromdichten und den Meßwerten an den Meßorten eine lineare Beziehung besteht, die durch eine Lead- Field-Matrix beschrieben ist, die nur von der relativen Lage der Zellen und der Meßorte zueinander abhängig ist,
b) aus der Lead-Field-Matrix wird eine nach Moore-Penrose verallgemeinerte inverse Lead-Field-Matrix gebildet,
c) bei Abwesenheit von elektrophysiologischen Aktivitäten wird an den Meßorten mit dem Vielkanalmeßsystem in einem ersten Meßintervall an mehreren Zeitpunkten Rauschen ge messen und von den Rauschmeßwerten eine Rauschmatrix gebildet, wobei jede Zeile der Matrix einem Meßort und jede Spalte der Matrix einem Meßzeitpunkt zugeordnet ist,
d) aus der Rauschmatrix wird eine Rauschkovarianzmatrix gebildet, die ein über das erste Meßintervall zeitlich gemitteltes Matrizenprodukt der Rauschmatrix mit der transponierten Rauschmatrix darstellt,
e) aus den innerhalb eines zweiten Meßintervalls an meh reren Zeitpunkten gemessenen Meßwerten des von elektro physiologischen Aktivitäten erzeugten Magnetfeldes wird eine Meßwertmatrix gebildet, wobei jede Zeile der Matrix einem Meßort und jede Spalte der Matrix einem Meßzeit punkt zugeordnet ist,
f) aus der Meßwertmatrix wird eine Meßwertkorrelations matrix gebildet, die ein über das zweite Meßintervall zeitlich gemitteltes Matrizenprodukt der Meßwertmatrix mit der transponierten Meßwertmatrix darstellt,
g) eine Differenzmatrix wird gebildet aus der Differenz der Meßwertekorrelationsmatrix und der Rauschkovarianzmatrix,
h) eine Stromdichtekorrelationsmatrix wird als Matrizen produkt aus der verallgemeinerten inversen Lead-Field- Matrix, der Differenzmatrix und der transponierten verallgemeinerten inversen Lead-Field-Matrix gebildet,
i) die Quadratwurzeln der diagonalen Elemente der Strom dichtekorellationsmatrix werden zur ortsrichtigen Dar stellung einer im zweiten Meßintervall gemittelten Stromdichteverteilung im Untersuchungsgebiet verwendet.

Eine Computer-Simulation hat gezeigt, daß selbst bei einem Signal-Rausch-Verhältnis von 0,2 noch eine Rekonstruktion der Stromdichteverteilung und damit eine Lokalisierung von elektrophysiologischen Aktivitäten möglich ist.

Ein Ausführungsbeispiel der Erfindung wird im folgenden anhand von 4 Figuren erläutert. Es zeigt

Fig. 1 den Aufbau eines biomagnetischen Meßsystems,

Fig. 2 ein Blockschaltbild des Verfahrens zur Lokalisie rung von elektrophysiologischen Aktivitäten, denen starkes Rauschen überlagert ist,

Fig. 3 eine zeitlich gemittelte Rekonstruktion eines ro tierenden und eines oszillierenden Dipols aus ver rauschten Meßwerten,

Fig. 4 zum Vergleich eine konventionelle Rekonstruktion der in Fig. 3 rekonstruierten Dipole aus verrauschten Meßwerten.

Zum Detektieren des von elektrophysiologischen Aktivitäten erzeugten sehr schwachen Magnetfeldes an Meßpunkten außer halb eines Untersuchungsgebietes wird ein biomagnetisches Meßsystem eingesetzt, dessen grundsätzlicher Aufbau in Fig. 1 angegeben ist. Fig. 1 zeigt schematisch eine magnetische Abschirmkammer 2, mit der außerhalb erzeugte Störfelder zum größten Teil abgeschirmt werden. Ein zu untersuchender Patient 4 befindet sich auf einer in der Abschirmkammer 2 angeordneten Patientenliege 6. Elektrophysiologische Akti vitäten, die hier durch einen Pfeil 8 symbolisiert sind, erzeugen ein elektrisches und magnetisches Feld, wobei hier nur das magnetische Feld 10 ausgewertet wird. Dazu wird das Magnetfeld mit einer Vielkanalmeßanordnung 12 oberhalb des Patienten 4 gemessen. Die Vielkanalmeßanordnung 12 umfaßt eine Vielkanal-Gradiometeranordnung 14 mit räumlich ge trennt angeordneten Gradiometern, die lediglich den Gra dienten der Magnetfeldverteilung erfassen und somit schon bei der Messung homogene Störfelder unterdrücken. Hier ist aus Gründen der Übersichtlichkeit eine Vielkanal-Gradio meteranordnung 14 mit fünfzehn einzelnen Gradiometern dargestellt, jedoch werden in der Praxis Vielkanal-Gradio meteranordnungen 14 mit mehr als 30 Kanälen, z. B. 37 Ka näle, verwendet. Die Gradiometer in Vielkanal-Gradiometer anordnung 14 sind jeweils mit einem SQUID (Super Conducting Quantum Interference Device) verbunden. Die Vielkanal-SQUID- Anordnung 16 und die Vielkanal-Gradiometeranordnung 14 sind in einem Kryostaten angeordnet und dort auf so niedrige Temperatur gehalten, daß Supraleitung vorherrscht.

Die Vielkanalmeßanordnung 12 ist mittels eines Stativs in einer Untersuchungsposition arretierbar. Die Untersuchungs position gibt die Meßorte der Gradiometer vor. In der dar gestellten Untersuchungsposition wird das Feld von zere bralen Aktivitäten gemessen. Die an den Meßorten zeitlich nacheinander gemessenen Meßsignale werden an eine Signal auswertungseinheit 18 übergeben, die sowohl das zeitliche Verhalten der Meßsignale anzeigt, als auch für ausgewählte Zeitpunkte eine äquivalente Stromdichteverteilung bestimmt, deren theoretisches Feld dem gemessenen Feld am nächsten kommt. Ein vollständiges Ersatzmodell besteht aus dem Ort, der Stärke und der Richtung der Stromdichte. Zum Ersatz modell gehört ebenfalls der Raum, in dem die Stromdichte verteilung angenommen wird. Der Raum, worin sich die äqui valente Stromdichteverteilung findet, ist in einem Ersatz modell für zerebrale Aktivitäten meist eine Kugel mit ho mogener Leitfähigkeit und für kardiologische Aktivitäten meist ein unendlicher Halbraum mit homogener Leitfähigkeit.

Bei der Lokalisierung der über ein Meßintervall zeitlich gemittelten elektrophysiologischen Aktivitäten über eine Rekonstruktion einer Stromdichteverteilung wird das Unter suchungsgebiet in N Zellen aufgeteilt, in denen jeweils eine Stromdichte j angenommen wird, die sich bei den Mo dellen "Kugel" und "Halbraum" aus maximal zwei Komponenten zusammensetzt. Die Gesamtzahl 2N der Komponenten ist größer als die Anzahl der Meßorte M. Zwischen den Stromdichten j in den Zellen und den Meßwerten des davon erzeugten mag netischen Feldes an den Meßorten besteht eine lineare Be ziehung, die durch eine Lead-Field-Matrix L beschrieben ist, die nur von der relativen Lage der Zellen und der Meßorte zueinander abhängig ist. Da 2N größer als M ist, handelt sich somit um ein unterbestimmtes System von Gleichungen zur Bestimmung der Stromdichten j aus dem an den Meßorten gemessenen Magnetfeld. Die Lead-Field-Matrix L ist durch die Meßposition x, y, z des Vielkanal-Meßsystems 12 vorgeben. Aus der Lead-Field-Matrix L wird eine nach Moore-Penrose verallgemeinerte inverse Lead-Field-Matrix L⁻ gebildet.

Bei Abwesenheit von elektrophysiologischen Aktivitäten wird an den Meßorten mit dem Vielkanalmeßsystem 12 in einem ersten Meßintervall an mehrere Zeitpunkten Rauschen ge messen und von den Rauschmeßwerten eine Rauschmatrix n ge bildet. Dabei ist jede Zeile der Matrix n einem Meßort m mit m=1 bis M und jede Spalte der Matrix n einem Meßzeit punkt t′_i mit i=1 bis K′ zugeordnet.

Aus der Rauschmatrix n wird eine Rauschkovarianzmatrix C gebildet, die ein über das Meßintervall zeitlich gemittel tes Matrizenprodukt der Rauschmatrix n mit der transpo nierten Rauschmatrix n^T darstellt. Die Rauschkovarianzma trix lautet demnach

C = <n n^T<.

Die spitze Klammer << bedeutet den zeitlichen Mittelwert.

Innerhalb eines zweiten Meßintervalls wird an mehreren Zeitpunkten das von elektrophysiologischen Aktivitäten er zeugte Magnetfeld gemessen und eine Meßwertmatrix B ge bildet. Dabei ist jede Zeile der Matrix B einem Meßort m mit m=1 bis M, und jede Spalte der Matrix B einem Meß zeitpunkt t_i mit i=1 bis K zugeordnet.

Aus der Meßwertmatrix B wird eine Meßwertkorrelationsma trix D als ein über das zweite Meßintervall zeitlich ge mitteltes Matrizenprodukt der Meßwertmatrix B mit der transponierten Meßwertmatrix B^T gebildet. Somit lautet die Meßwertkorrelationsmatrix

D = <B B^T<.

Die spitzen Klammern << symbolisieren die zeitliche Mittelwertbildung.

Eine Differenzmatrix Z wird gebildet aus der Differenz der Meßwertkorrelationsmatrix D und der Rauschkovarianzmatrix C. Somit ist die Differenzmatrix

Z = D-C.

Eine Stromdichtekorrelationsmatrix S wird als Matrizen produkt aus der verallgemeinerten inversen Lead-Field- Matrix L⁻, der Differenzmatrix Z und der transponierten verallgemeinerten inversen Lead-Field-Matrix (L^-1)^T ge bildet.

Die Quadratwurzeln der diagonalen Elemente der Stromdichte korrelationsmatrix S werden zur ortsrichtigen Darstellung einer im zweiten Meßintervall gemittelten Stromdichtever teilung im Untersuchungsgebiet verwendet.

Dazu werden die Stromdichten nach Größe und Richtung sowie entsprechend den Koordinaten der Zellen auf einem Bild schirm 20 dargestellt. Als Beispiel für die ortsrichtige Darstellung wird auf die weiter unten beschriebenen Fig. 3 und 4 verwiesen.

Zur Überprüfung der Leistungsfähigkeit des Lokalisierungs verfahrens wurden Computersimulationen durchgeführt, wobei den Meßwerten Rauschen überlagert wurde.

Bei der Erzeugung der Meßwertmatrix wurde von zwei einzel nen Dipolen ausgegangen, die 5 cm unterhalb des Vielkanal meßsystems 12 lokalisiert waren. Dabei wurde ein unend licher leitfähiger Halbraum angenommen. Der erste Dipol rotierte in der x-y-Ebene am Ort (in kartesischen Koordi naten) x=-3 cm, y=-1 cm mit einer Stärke von 1 mA mm. Der zweite Dipol war in derselben Ebene wie der erste Dipol angeordnet mit den Koordinaten x=4 cm , y=2,5 cm. Das Dipolmoment in x-Richtung war Null, während das Dipolmoment in y-Richtung sinusförmig moduliert wurde. Die Amplitude des Dipolmomentes war Dy=1,2 mA mm.

Das Verhältnis der gemittelten Signalleistung

auf die Rauschleistung bezogen, also

war 0,2. Die Standardabweichung des Rauschens war 16 pT.

Bei der Rekonstruktion wurde die Meßwertkorrelationsmatrix D aus den verrauschten Meßwerten B für 20 000 Meßzeitpunkte bestimmt. Die Rekonstruktion der Stromdichteverteilung j wurde in einer Ebene bestimmt, die exakt der Tiefe der Stromdipole entsprach. Die Rekonstruktion wurde entspre chend wie anhand von Fig. 2 beschrieben durchgeführt. Das Ergebnis der Rekonstruktion zeigt Fig. 3. Die Stromdichten j liegen ausnahmslos in der x-y-Ebene und haben nur eine x- und y-Komponente jx bzw. jy. Die Pfeile geben die Rich tungen der Stromdichten an den Orten oder Zellen des Un tersuchungsgebiets an, an denen sich die Pfeile befinden. Die Größe der Pfeile ist ein Maß für den Betrag der Strom dichte an dem Ort.

Zum Vergleich wurde eine konventionelle Rekonstruktion durchgeführt, bei der der Einfluß des Rauschens nicht durch die Einführung der Rauschkovarianzmatrix C und der Meßwert korrelationsmatrix D in das Lokalisierungsverfahrens ver ringert wurde. Die konventionelle Rekonstruktion kann nur für die Meßwerte eines einzigen Zeitpunktes durchgeführt werden. Die zugrundeliegenden Werte des ersten Dipols waren Dx=1/√mA mm und Dy=1/√mA mm. Der zweite Dipol hatte ein konstantes Moment von Dy=1 mA mm. Diesen Werten wurde Rauschen überlagert. Die Ergebnisse der konventio nellen Rekonstruktion zeigt Fig. 4.

Vergleicht man die Ergebnisse, die durch das konventionelle Rekonstruktionsverfahren gewonnen wurden, mit den Ergeb nissen des zeitgemittelten Rekonstruktionsverfahrens unter Berücksichtigung der Rauschkovarianzmatrix C, so erkennt man, daß durch die Verwendung der Rauschkovarianzmatrix C und der Meßwertkorrelationsmatrix D der negative Einfluß des externen Rauschens auf die Genauigkeit fast vollständig beseitigt wurde. Dagegen zeigt die Fig. 4, daß bei stark verrauschten Meßwerten die konventionelle Rekonstruktion nahezu unbrauchbare Ergebnisse liefert.

Claims

Verfahren zur Lokalisierung von elektrophysiologischen Aktivitäten (8), die in einem Lebewesen (4) innerhalb eines Untersuchungsgebiets auftreten und die ein Magnetfeld (10) erzeugen, das an Meßpunkten außerhalb des Untersuchungs gebiets mit einem Vielkanalmeßsystem (12) gemessen wird, mit den Schritten:
- a) in einem Modell der elektrophysiologischen Aktivitäten im Untersuchungsgebiet wird das Untersuchungsgebiet in Zellen aufgeteilt, in denen jeweils eine Stromdichte (j) vorhanden ist, die sich aus maximal drei Komponenten zusammensetzt, wobei die Gesamtzahl der Komponenten der größer ist als die Anzahl der Meßpunkte (M) und wobei zwischen den Stromdichten (j) und den Meßwerten (B) an den Meßorten eine lineare Beziehung besteht, die durch eine Lead-Field-Matrix (L) beschrieben ist, die nur von der relativen Lage der Zellen und der Meßorte zueinander abhängig ist,
- b) aus der Lead-Field-Matrix wird eine nach Moore-Penrose verallgemeinerte inverse Lead-Field-Matrix (L⁻) gebil det,
- c) bei Abwesenheit von elektrophysiologischen Aktivitäten wird an den Meßorten mit dem Vielkanalmeßsystem (12) in einem ersten Meßintervall an mehreren Zeitpunkten Rauschen gemessen und von den Rauschmeßwerten eine Rauschmatrix (n) gebildet, wobei jede Zeile der Matrix einem Meßort und jede Spalte der Matrix einem Meßzeit punkt (t′_i) zugeordnet ist,
- d) aus der Rauschmatrix wird eine Rauschkovarianzmatrix gebildet, die ein über das erste Meßintervall zeitlich gemitteltes Matrizenprodukt der Rauschmatrix (n) mit der transponierten Rauschmatrix (n^T) darstellt,
- e) aus den innerhalb eines zweiten Meßintervalls an mehreren Zeitpunkten (t_i) gemessenen Meßwerten des von elektrophysiologischen Aktivitäten erzeugten Magnet feldes wird eine Meßwertmatrix (B) gebildet, wobei jede Zeile der Matrix einem Meßort und jede Spalte der Matrix einem Meßzeitpunkt (t_i) zugeordnet ist,
- f) aus der Meßwertmatrix wird eine Meßwertkorrelations matrix (D) gebildet, die ein über das zweite Meßinter vall zeitlich gemitteltes Matrizenprodukt der Meßwert matrix (B) mit der transponierten Meßwertmatrix (B^T) darstellt,
- g) eine Differenzmatrix (Z) wird gebildet aus der Differenz der Meßwertkorrelationsmatrix (D) und der Rauschko varianzmatrix (C)
- h) eine Stromdichtekorrelationsmatrix (S) wird als Matri zenprodukt aus der verallgemeinerten inversen Lead- Field-Matrix (L⁻), der Differenzmatrix und der trans ponierten verallgemeinerten inversen Lead-Field-Matrix ((L⁻)^T) gebildet
- i) die Quadratwurzeln der diagonalen Elemente der Strom dichtekorrelationsmatrix (S) werden zur ortsrichtigen Darstellung einer im zweiten Meßintervall gemittelten Stromdichteverteilung im Untersuchungsgebiet verwendet.