KR101525485B1 - 3차원 심장 윤곽 재구성 방법 - Google Patents

3차원 심장 윤곽 재구성 방법 Download PDF

Info

Publication number
KR101525485B1
KR101525485B1 KR1020140033636A KR20140033636A KR101525485B1 KR 101525485 B1 KR101525485 B1 KR 101525485B1 KR 1020140033636 A KR1020140033636 A KR 1020140033636A KR 20140033636 A KR20140033636 A KR 20140033636A KR 101525485 B1 KR101525485 B1 KR 101525485B1
Authority
KR
South Korea
Prior art keywords
heart
signal
source current
contour
matrix
Prior art date
Application number
KR1020140033636A
Other languages
English (en)
Inventor
김기웅
하태훈
권혁찬
이용호
김진목
Original Assignee
한국표준과학연구원
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 한국표준과학연구원 filed Critical 한국표준과학연구원
Priority to KR1020140033636A priority Critical patent/KR101525485B1/ko
Priority to PCT/KR2015/002569 priority patent/WO2015142029A1/ko
Priority to CN201580015416.3A priority patent/CN106132288B/zh
Application granted granted Critical
Publication of KR101525485B1 publication Critical patent/KR101525485B1/ko

Links

Images

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/24Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
    • A61B5/242Detecting biomagnetic fields, e.g. magnetic fields produced by bioelectric currents
    • A61B5/243Detecting biomagnetic fields, e.g. magnetic fields produced by bioelectric currents specially adapted for magnetocardiographic [MCG] signals

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Biophysics (AREA)
  • Pathology (AREA)
  • Engineering & Computer Science (AREA)
  • Biomedical Technology (AREA)
  • Cardiology (AREA)
  • Physics & Mathematics (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • General Health & Medical Sciences (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

본 발명은 심장 윤곽 재구성 방법 및 심자도 측정 시스템을 제공한다. 이 심장 윤곽 재구성 방법은 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 제1 윤곽을 구하는 단계; 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계; 및 상기 제1 윤곽과 상기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성하는 단계를 포함한다.

Description

3차원 심장 윤곽 재구성 방법{Three-dimensional Reconstruction of a Cardiac Outline by Magnetocardiography}
본 발명은 심장 윤곽 재구성 방법에 관한 것으로, MRI 장치와 같은 해부학적 영상 장치의 도움 없이 심자도 장치로 측정한 다채널 데이터를 이용하여 심장 윤곽을 3차원으로 형성하는 방법에 관한 것이다.
심자도(Magnetocardiography; MCG)는 초고감도 자기센서인 초전도양자간섭장치 (Superconducting Quantum Interference Device; SQUID) 센서를 통해 심근(myocardium)의 전기활동(electric activity)에 의해 발생되는 자기 신호를 측정하여 심장질환을 비침습적으로 진단하는 장치이다. 심자도 장치로 측정한 다채널 데이터를 이용하면 심장질환의 병소를 국지화 할 수 있다. 특히, 3차원 심장모형에 심장질환의 병소를 국지화하면 의사에게 병소의 정보를 직관적으로 보여줄 수 있다.
하지만 상기 3차원 심장모형을 얻기 위해서는 환자의 엑스선 컴퓨터단층촬영 (X-ray computed tomography; CT)장치 또는 자기공명영상(magnetic resonance imaging; MRI)장치로 측정한 해부학적 영상데이터가 필요하다. 상기 해부학적 영상데이터를 얻기 위하여, 환자는 불필요한 방사선 피폭이나 강한 자기장에 노출될 수 있다. 또한, 심장모형을 얻기 위해서는 해부학적 영상데이터로 심장을 분리(segmentation)해내는 과정이 필요하다. 상기 심장분리과정은 전문가가 수동으로 작업하거나 또는 자동화된 심장분리 알고리즘을 사용할 수 있지만, 많은 시간이 소요된다.
Kenji Nakai 등에 의한 미국공개특허 제2008-0033312호는 비적응 공간 필터링(non-adaptive spatial filtering) 방법 중 하나인 최소크기추정법(minimum norm estimation)에 티코노프 정규화(Tikhonov regularization) 방법을 적용하여 전류 밀도를 추정한 후 이를 바탕으로 3차원 심장 윤곽을 생성하는 방법을 보고하였다. 하지만, 비적응 공간 필터링 방법은 깊이 있는 신호원에 대한 추정이 정확하지 않기 때문에, 3차원 심장 모형을 정확하게 생성하는데 한계가 있다.
본 발명의 해결하고자 하는 일 기술적 과제는 3차원 심장 시각화 또는 매핑 방법(3-D cardiac visualization or mapping method)을 제공하는 것이다.
본 발명의 일 실시예에 따른 심장 윤곽 재구성 방법은 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 제1 윤곽을 구하는 단계; 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계; 및 상기 제1 윤곽과 상기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성하는 단계를 포함한다.
본 발명의 일 실시예에 있어서, 심장의 제1 윤곽을 구하는 단계는: 심자도 장치를 이용하여 심근에서 발생된 자기 신호를 측정하는 단계; 다채널로 측정된 심자도 신호에서 심방과 심실을 나타내는 P파와 T파 구간의 심장 자기 신호를 추출하는 단계; 심장의 단위 소스전류 및 심장을 근사화한 평면체적도체모델 (horizontally layered volume conductor model)과 심자도 센서의 위치정보를 이용하여 자기장(리드필드 벡터)를 계산하는 단계; 및 거리에 따른 신호원의 크기를 보상하기 위하여 상기 리드필드 벡터를 규격화(normalization)하는 단계를 포함할 수 있다.
본 발명의 일 실시예에 있어서, 심장의 제1 윤곽을 구하는 단계는: 임의의 가중치행렬을 단위 행렬로 초기화하는 단계; 상기 리드필드 벡터와 임의의 가중치행렬을 조합하여 그램행렬을 구하는 단계; 상기 구해진 그램행렬에 측정신호의 신호대 잡음비를 고려하여 정규화 매개변수(regularization parameter)를 설정하는 단계; 상기 정규화된 그램행렬과 규격화된 리드필드 벡터를 조합하여 배열이득제한최소크기 공간필터(array-gain constraint minimum-norm spatial filter)를 구하는 단계; 상기 구해진 배열이득제한최소크기 공간필터와 측정된 심방과 심실 신호를 이용하여 추정 소스전류 벡터의 파워를 구하는 단계; 상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하는지 확인하는 단계; 상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하지 않는 경우 가중치행렬을 추정 소스벡터를 이용하여 갱신하는 단계; 상기 구해진 추정 소스 전류 벡터의 파워를 규격화하는 단계; 및 규격화된 소스전류벡터의 파워를 신호대잡음비로 설정된 제1 임계값을 적용하여 심방과 심실의 제1 윤곽을 재구성하는 단계들 중에서 적어도 하나를 더 포함할 수 있다.
본 발명의 일 실시예에 있어서, 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계는: 심방과 심실에서 소스전류벡터의 파워가 가장 큰 지점을 각각 구하는 단계; 측정된 심자도 신호와 회귀갱신그램행렬 배열이득제한최소크기 공간필터를 선형결합하여 모든 소스영역에서 신호파형(추정 소스전류 벡터)을 재구성하는 단계; 상기 재구성된 신호파형(추정 소스전류 벡터)을 고속 퓨리에 변환(fast Fourier transform)하는 단계; 상기 고속 퓨리에 변환된 신호를 이용하여 심방과 심실에서 모든 소스영역에 대한 상호-스펙트럼(cross-spectrum)과 자기-스펙트럼(auto-spectrum)을 구하는 단계; 상기 구해진 상호-스펙트럼과 자기-스펙트럼을 조합하여 심방의 결맞음과 심실의 결맞음을 구하는 단계; 상기 구해진 심방과 심실의 결맞음을 규격화하는 단계; 및 상기 규격화된 결맞음을 신호대 잡음비로 설정된 제2 임계값을 적용하여 심방과 심실의 제2 윤곽을 재구성하는 단계들 중에서 적어도 하나를 포함할 수 있다.
본 발명의 일 실시예에 따른 심장 윤곽 재구성 방법은 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 윤곽을 구한다.
본 발명의 일 실시예에 있어서, 심장의 윤곽을 구하는 단계는 심자도 장치를 이용하여 심근에서 발생된 자기 신호를 측정하는 단계; 다채널로 측정된 심자도 신호에서 심방과 심실을 나타내는 P파와 T파 구간의 심장 자기 신호를 추출하는 단계; 심장의 단위 소스전류 및 심장을 근사화한 평면체적도체모델 (horizontally layered volume conductor model)과 심자도 센서의 위치정보를 이용하여 자기장(리드필드 벡터)를 계산하는 단계; 및 거리에 따른 신호원의 크기를 보상하기 위하여 상기 리드필드 벡터를 규격화(normalization)하는 단계를 포함할 수 있다.
본 발명의 일 실시예에 있어서, 임의의 가중치행렬을 단위 행렬로 초기화하는 단계; 상기 리드필드 벡터와 임의의 가중치행렬을 조합하여 그램행렬을 구하는 단계; 상기 구해진 그램행렬에 측정신호의 신호대 잡음비를 고려하여 정규화 매개변수(regularization parameter)를 추가하는 단계; 상기 정규화된 그램행렬과 규격화된 리드필드 벡터를 조합하여 배열이득제한최소크기 공간필터(array-gain constraint minimum-norm spatial filter)를 구하는 단계; 상기 구해진 배열이득제한최소크기 공간필터와 측정된 심방과 심실 신호를 이용하여 추정 소스전류 벡터의 파워를 구하는 단계; 상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하는지 확인하는 단계; 상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하지 않는 경우 가중치행렬을 추정 소스벡터를 이용하여 갱신하는 단계; 상기 구해진 추정 소스 전류 벡터의 파워를 규격화하는 단계; 규격화된 소스전류벡터의 파워를 신호대잡음비로 설정된 임계값을 적용하여 심방과 심실의 윤곽을 재구성하는 단계를 더 포함할 수 있다.
본 발명의 일 실시예에 따른 심장 윤곽 재구성 방법은 심방과 심실에서 소스전류벡터의 파워가 가장 큰 지점을 각각 구하는 단계; 측정된 심자도 신호와 회귀갱신그램행렬 배열이득제한최소크기 공간필터를 선형결합하여 모든 소스영역에서 신호파형(추정 소스전류 벡터)을 재구성하는 단계; 상기 재구성된 신호파형(추정 소스전류 벡터)을 고속 퓨리에 변환(fast Fourier transform)하는 단계; 상기 고속 퓨리에 변환된 신호를 이용하여 심방과 심실에서 모든 소스영역에 대한 상호-스펙트럼(cross-spectrum)과 자기-스펙트럼(auto-spectrum)을 구하는 단계; 상기 구해진 상호-스펙트럼과 자기-스펙트럼을 조합하여 심방의 결맞음과 심실의 결맞음을 구하는 단계; 상기 구해진 심방과 심실의 결맞음을 규격화하는 단계; 및 상기 규격화된 결맞음을 신호대 잡음비로 설정된 임계값을 적용하여 심방과 심실의 윤곽을 재구성하는 단계를 포함한다.
본 발명의 일 실시예에 따른 프로그램을 기록한 기록매체는 위의 방법들을 컴퓨터에서 실행시킬 수 있다.
본 발명의 일 실시예에 따른 심자도 측정 시스템은 심자도 센서 및 자기차폐실을 포함하고 심자도 신호를 측정하는 심자도 측정 장치; 및 상기 심자도 신호를 처리하여 심장의 윤곽을 재구성하는 처리부를 포함한다. 상기 처리부는 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 제1 윤곽을 구한다. 상기 처리부는 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구한다. 상기 처리부는 상기 제1 윤곽과 상기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성한다.
본 발명의 일 실시예에 따른 심장 윤곽 3차원 재구성 방법은, 해부학적 영상을 얻기 위해 필요한 CT나 MRI 장치의 사용 없이 환자의 MCG 데이터만을 이용해 심장모형을 재구성하는 것이 가능하다. 따라서 환자는 불필요한 방사선 피폭이나 강한 자기장에 노출되는 위험으로부터 피할 수 있다. 또한, 상기 영상장치를 사용하지 않을 경우 환자의 경제적인 부담도 줄일 수 있다. 추가적으로, 심장모형을 만들기 위한 심장분리과정이 필요하지 않기 때문에 MCG 장치를 사용하여 실시간으로 환자를 진단하는데 활용될 수 있다.
도 1은 본 발명의 일 실시예에 따른 3차원 심장 윤곽 재구성 방법을 설명하는 흐름도이다.
도 2a, 도 2b, 및 도 2c는 본 발명의 일 실시예에 따른 3차원 심장 윤곽 재구성 방법을 설명하는 흐름도이다.
도 3은 본 발명의 일 실시예에 따른 심자도 측정 장치를 설명하는 도면이다.
도 4는 본 발명의 일 실시예에 따른 다채널로 측정된 심자도 신호에서 심장파형 추출을 나타내는 도면이다.
도 5는 본 발명의 일 실시예에 따른 심자도의 리드필드벡터(lead-field vector) 계산에 필요한 소스 공간과 센서 공간의 관계를 나타내는 도면이다.
도 6은 본 발명의 일 실시예에 따른 심방의 결맞음 매핑을 나타내는 도면이다.
도 7은 본 발명의 일 실시예에 따른 심실의 결맞음 매핑을 나타내는 도면이다.
도 8은 본 발명의 일 실시예에 따른 수치 시뮬레이션을 위한 가상 소스 공간 및 가상 심장을 나타내는 도면이다.
도 9는 본 발명의 일 실시예에 따른 수치 시뮬레이션의 결과를 나타내는 도면이다.
도 10은 본 발명의 일 실시예에 따른 모형 실험을 위한 모형 심장(cardiac phantom)을 나타내는 도면이다.
도 11은 도 10의 모형 심장을 이용한 수치 시뮬레이션 결과를 나타내는 도면이다.
도 12는 신호대 잡음비에 따른 정규화 매개변수의 비율을 나타내는 도면이다.
3차원 심장 시각화 또는 매핑 방법(3-D cardiac visualization or mapping method)은 심자도(magnetocardiography;MCG) 임상적 응용(clinical applications)에 도움이 된다. 그러나, 심장 재구성(cardiac reconstruction)은 추가적인 영상 양식(image modalities)을 요구한다. 우리는 추가적인 영상 기술 없이 MCG 측정 데이터만을 사용하여 3차원 심장 윤곽 재구성 방법을 제안한다. 심장 윤곽은 공간 필터링(spatial filtering) 방법과 결맞음 매핑(coherence mapping) 방법의 조합에 의하여 재구성될 수 있다.
심장 활동성의 세기(The strength of cardiac activities)는 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간 필터에 의하여 표시될 수 있다.
또한, 심방의 최대 소스 점과 모든 소스 점 사이의 결맞음과 심실의 최대 소소점과 모든 소스 점 사이의 결맞음은 결맞음 매빙 방법에 의하여 비교되었다.
수치 시뮬레이션과 인체 모형 실험(phantom experiments)은 원래의 모양과 재구성된 구성을 비교하여 3차원 심장 윤곽 재구성의 효율성을 증명하였다.
본 발명의 일 실시예에 따른 심장 윤곽 3차원 재구성 방법은 CT 또는 MRI 장치로부터 얻은 해부학적 영상을 사용하지 않고, 환자로부터 측정된 심자도 (MCG) 데이터만을 사용하여 심장 윤곽을 3차원으로 생성하는 방법이다. 구체적으로 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류 파워(source current power)로부터 심장의 제1 윤곽이 구해진다. 또한, 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽이 구해진다. 마지막으로 제1 윤곽과 제2 윤곽을 상호 결합하여 심장의 제3 윤곽이 구성된다.
심자도(MCG)는 심근에서 발생된 자기장을 비침습적으로 측정할 수 있는 장치이다. 심자도는 시간 및 공간 분해능(time and spatial resolution)이 뛰어나기 때문에 심근허혈 영역(myocardial ischemic regions) 또는 심부정맥(cardiac arrhythmia)의 원인이 되는 회귀성 흥분 영역(reentrant excitation regions)을 추정하는데 아주 유용하다. 최근 심자도를 이용하여 심근의 회귀성 흥분 영역을 3차원 심장 표면에 가시화시키는 연구가 보고되었다. 심장 전기 활동(cardiac electric activities)을 3차원 심장 모형에 가시화하면 임상적으로 많은 응용이 가능하다.
그러나 통상적으로 3차원 심장 모형을 얻기 위해서는 몇 가지 단점이 존재한다. 첫째, 엑스선 컴퓨터 단층 촬영(CT) 또는 자기공명영상(MRI) 등과 같은 해부학적 영상 장치들로부터 복잡한 분할(segmentation) 과정이 요구된다. 최근 다수의 연구 그룹이 수동 분할과정의 단점을 극복한 자동 분할 방법을 제안하였다. 하지만, 여전히 해부학적 영상 취득과정이 반드시 필요하다. 둘째로, 해부학적 영상정보를 얻을 때 환자는 엑스선이나 높은 자기장에 노출되며, 이러한 노출은 잠재적인 위험을 초래할 수 있다. 마지막으로 해부학적 영상정보를 얻기 위해 환자는 비싼 비용을 지불해야 된다.
상기 3차원 심장 모형 획득과정의 단점을 극복하기 위하여, 심자도 측정 데이터만을 사용하여 3차원 심장 모형을 생성하는 방법이 제안되었다. 일본의 나카기(Nakai) 그룹은 비적응 공간 필터링(non-adaptive spatial filtering) 방법 중 하나인 최소크기추정법(minimum norm estimation)에 티코노프 정규화(Tikhonov regularization) 방법을 적용하여 전류 밀도를 추정한 후 이를 바탕으로 3차원 심장 윤곽을 생성하는 방법을 보고하였다. 하지만, 비적응 공간 필터링 방법은 깊이 있는 신호원에 대한 추정이 정확하지 않는다. 따라서, 비적응 공간 필터링 방법은 3차원 심장 모형을 정확하게 생성하는데 한계가 있다.
비적응 공간필터링에 비하여, 적응(adaptive) 공간필터링은 높은 공간분해능을 보여준다. 적응 공간필터링은 측정대상의 기하학적 형태(measurement geometry)와 공분산행렬(covariance matrix)를 이용하기 때문에 노이즈간섭(noise interferences)에 강하다. 공간필터링은 약하게 상관된(correlated) 신호에는 강하지만, 강하게 상관된 신호에 대해서는 큰 국지화(localization) 오차를 나타낸다. 심장의 전기신호는 강하게 상관되어있다. 따라서, 적응 공간필터링으로 심근전류원을 국지화하면 큰 오차가 나타난다. 따라서 비적응 공간필터링과 적응 공간필터링 모두 3차원 심장 모형을 재구성하는데 한계점이 있다.
최근 공간해상도를 높인 회귀갱신그램행렬 배열이득제한최소크기 공간필터가 제안되었다. 상기 공간필터는 비적응 공간필터를 기반으로 적응 공간필터에 상응하는 결과를 도출할 수 있다. 따라서 상기 공간필터를 사용하면 강하게 상관된 심근 신호에서 발생된 문제를 회피함으로써 더 정확하게 심장 모형을 생성할 수 있다.
또한 우리는 심장 모형을 더 정확하게 생성하기 위해 결맞음 매핑 방법을 추가적으로 적용하였다. 심방과 심실에서 생성된 전기신호는 각각의 영역에서 서로 상관된 파형을 보여준다. 따라서 재생성된 심장파형의 비슷함을 비교함으로써 심방과 심실의 더 정확한 모형을 재구성할 수 있다.
이하, 첨부된 도면을 참조로 본 발명의 바람직한 실시예들에 대하여 보다 상세히 설명한다. 이하의 도면들에서 동일한 참조부호는 동일한 구성요소를 지칭하며, 도면상에서 각 구성요소의 크기는 설명의 명료성과 편의상 과장되어 있을 수 있다.
도 1은 본 발명의 일 실시예에 따른 3차원 심장 윤곽 재구성 방법을 설명하는 흐름도이다.
도 2a, 도 2b, 및 도 2c는 본 발명의 일 실시예에 따른 3차원 심장 윤곽 재구성 방법을 설명하는 흐름도이다.
도 1 및 도 2a 내지 도 2c를 참조하면, 3차원 심장 윤곽 재구성 방법은 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류 파워(source current power)로부터 심장의 제1 윤곽을 구하는 단계(S100); 결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계(S200); 및 상기 제1 윤곽과 싱기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성하는 단계(S300)를 포함한다.
심장의 제1 윤곽을 구하는 단계(S100)는 다음과 같이 처리된다. 심자도 장치를 이용하여 심근에서 발생된 심자도 신호가 측정된다(S110). 다채널로 측정된 심자도 신호에서 심방과 심실을 나타내는 P파와 T파 구간의 심장 자기 신호(ba(t),bv(t))가 추출된다(S120). 심장의 단위 소스전류 및 심장을 근사화한 평면체적도체모델 (horizontally layered volume conductor model)과 심자도 센서의 위치정보를 이용하여 자기장(리드필드 벡터;L(r))가 계산된다(S130). 거리에 따른 신호원의 크기를 보상하기 위하여 상기 리드필드 벡터가 규격화(normalization)된다(S140).
회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터가 구해지는 방법이 소개된다. 임의의 가중치행렬(H(r))는 단위 행렬(I)로 초기화된다(S150). 상기 리드필드 벡터(L(r))와 임의의 가중치행렬(H(r))을 조합하여 그램행렬(G)이 구해진다(S160). 상기 구해진 그램행렬(G)에 측정신호의 신호대 잡음비를 고려하여 정규화 매개변수(regularization parameter;λ)가 설정된다(S170).
정규화(regulariztion)된 그램행렬(G)과 규격화된 리드필드 벡터를 조합하여 배열이득제한최소크기 공간필터(array-gain constraint minimum-norm spatial filter)가 구해진다(S180). 상기 구해진 배열이득제한최소크기 공간필터와 측정된 심방과 심실 신호를 이용하여 추정 소스전류 벡터 및 추정 소스 전류의 파워가 구해진다(S191). 상기 구해진 추정 소스전류의 파워가 일정한 값으로 수렴하는지 확인된다(S192). 상기 구해진 추정 소스전류의 파워가 일정한 값으로 수렴하지 않는 경우 가중치 행렬는 상기 추정 소스전류 벡터를 이용하여 갱신된다(S193). 상기 구해진 추정 소스전류의 파워가 일정한 값으로 수렴하는 경우 상기 구해진 추정 소스 전류 벡터의 파워는 규격화(normalization)된다(S194). 규격화된 추정 소스전류의 파워를 신호대잡음비로 설정된 제1-1 임계값과 제1-2 임계값을 적용하여 심방과 심실의 제1 윤곽은 재구성된다(195).
결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계는 다음과 같다.
위에서 설명한 방법 또는 다른 방법으로, 심장의 추정 소스전류 벡터 또는 심장의 추정 소스전류의 파워가 계산된다. 심방과 심실에서 소스전류의 파워가 가장 큰 지점(ra, rb)이 각각 구해진다(S210). 측정된 심장 자기 신호와 회귀갱신그램행렬 배열이득제한최소크기 공간필터를 선형결합하여 모든 소스영역에서 소스전류신호파형(
Figure 112014027564777-pat00001
)이 재구성된다(S220). 본 발명의 변형된 실시예에 따르면, 상기 소스전류신호파형은 S191에서 계산된 값을 사용할 수 있다.
상기 재구성된 소스전류 신호파형은 고속 퓨리에 변환(fast Fourier transform)된다(S230). 상기 고속 퓨리에 변환된 신호를 이용하여 심방과 심실에서 모든 소스영역에 대한 상호-스펙트럼(cross-spectrum; Γaq(f),Γvq(f))과 자기-스펙트럼(auto-spectrum;Γaa(f),Γvv(f))이 구해진다(S240). 상기 심방 상호-스펙트럼은 심방에 소스전류의 파워가 가장 큰 지점과 다른 지점 사이의 상호 상관(corss-correlation)이다. 상기 심실 상호-스펙트럼은 심실에 소스전류의 파워가 가장 큰 지점과 다른 지점 사이의 상호 상관(corss-correlation)이다.
상기 심방 자기-스펙트럼은 심방에 소스전류의 파워가 가장 큰 지점의 자기-상관(auto-correlation)이다. 상기 심실 자기-스펙트럼은 심실에 소스전류의 파워가 가장 큰 지점의 자기-상관(auto-correlation)이다.
상기 구해진 상호-스펙트럼과 자기-스펙트럼을 조합하여 심방의 결맞음과 심실의 결맞음(Cohaq(f),Cohvq(f))이 구해진다(S250). 상기 결맞음은 규격화될 수 있다. 상기 규격화된 결맞음을 신호대 잡음비로 설정된 제2-1 임계값과 제2-2 임계값을 적용하여 심방과 심실의 제2 윤곽은 재구성된다(S260).
이어서, 상기 제1 윤곽과 상기 제2 윤곽은 상호 결합되어 심장의 제3 윤곽을 구성한다(S300).
도 3은 본 발명의 일 실시예에 따른 심자도 측정 장치를 설명하는 도면이다.
도 3을 참조하면, 심자도 측정 장치(100)는 자기차폐실(magnetically shielded room;101)내에 설치된다. 상기 심자도 장치(100)는 64채널의 초전도양자간섭소자(Superconducting Quantum Interference Device, SQUID)를 포함한다. 픽업 코일(pickup coil)은 70 mm의 베이스라인(baseline)을 가지는 1차 권선형 미분계(1st order axial gradiometer)이다. 각 SQUID 센서(112)의 노이즈 레벨은 100 Hz에서 5 fTrms/Hz1/2이고 샘플링률(sampling rate)은 500 Hz이다. 상기 SQUID 센서(112)는 인체에서 발생하는 미세한 자기신호를 측정한다. 상기 SQUID 센서(112)는 영하 250도 이하의 극저온 상태에서 작동할 수 있다. 따라서 상기 SQUID 센서는 극저온을 유지시키는 냉각 수단(114) 내부에 설치된다. 상기 냉각 수단(114)은 듀어(dewar,116)내부에 배치될 수 있다.
구동회로(124)는 상기 심자도 장치를 구동한다. 증폭기 및 필터(118)는 자기자폐실 내부에 배치된다. 전원부는 증폭기 및 필터(118)에 등에 전력을 공급할 수 있다. 측정된 심자도 신호는 처리부(122)로 전송되어 신호처리 된 후 분석된다(S110).
도 4는 본 발명의 일 실시예에 따른 다채널로 측정된 심자도 신호에서 심장파형 추출을 나타내는 도면이다.
도 4를 참조하면, 측정된 심자도 신호에서 상기 심장파형은 P, Q, R, S, T 피크가 순차적으로 나타난다. P 피크를 포함하는 구간은 P 파로 분리되고, Q, R, S 피크를 포함하는 구간은 QRS 파로 분리되고, T 피크를 포함하는 구간은 T 파로 분리된다. 상기 심장파형 중 P 파(ba(t)) 신호는 심방에서 발생된 신호이고, QRS 파는 심실이 수축할 때 발생된 신호이고, T 파(bv(t))는 심실이 이완할 때 발생된 신호이다.
본 발명에서 심방모형을 생성하기 위해서 P 파의 시작과 끝 부분을 추출하고, 심실모형을 생성하기 위해서 T 파의 시작과 끝 부분을 추출한다(S120). 추출된 P 파(ba(t))와 T 파(bv(t))는 다음과 같이 주어진다.
Figure 112014027564777-pat00002
여기서 n은 심자도 센서의 개수이고, t는 심장 파형의 추출 구간이다. 윗 첨자 T는 행렬의 전치(transpose)이다.
도 5는 본 발명의 일 실시예에 따른 심자도의 리드필드벡터(lead-field vector) 계산에 필요한 소스 공간과 센서 공간의 관계를 나타내는 도면이다.
도 5를 참조하면, 리드필드벡터는 단위 소스전류(unit current source)와 도체모델(conductor model)을 알고 있을 때, 센서 공간에서 계산된 자기장이다. 상기 소스전류는 크기와 방향을 가지는 등가전류쌍극자(equivalent current dipole; ECD)가 사용된다. 인체는 일정한 전기전도도를 가지는 도체 공간(z<0)으로 간주될 수 있다. 상기 도체모델은 평면도체모델(horizontally layered conductor model)이 사용된다. 평면도체모델은 평면 아래(z<0)는 도전체로 가정하고 평면 위(z>0)는 절연체로 가정한다. 상기 소스전류는 도체 공간 내에 존재한다. 또한 상기 도체 공간 중에서 소스 공간(V)이 선택된다. 심장은 상기 소스 공간(V) 내에 존재한다.
상기 소스 공간은 xyz 직각 좌표계를 기준으로 직육면체일 수 있다. 상기 소스 공간은 10 mm 단위로 분할된 복셀(voxels)로 구성될 수 있다. 상기 소스 공간은 x축 방향으로 150 mm이고, y축 방향으로 170 mm이고, z축 방향으로 150 mm일 수 있다. 이에 따라, 총 복셀의 개수(Q)는 3314개 일 수 있다.
상기 도체 공간의 외부에는 측정 영역(z>0)이 배치된다. 상기 측정 영역에 복수의 심자도 센서(112)가 배치되는 센서 공간이 배치될 수 있다. 상기 심자도 센서(112)는 SQUID 센서일 수 있다. 상기 심자도 센서의 위치 정보는 이미 알고 있다. 상기 SQUID 센서(112)는 70 mm의 베이스라인(baseline)을 가지는 1차 권선형 미분계(1st order axial gradiometer)를 포함할 수 있다.
공간 필터의 출력 파워는 심근 전기 활성(mycocardial electric activity)의 크기를 추정할 수 있기 때문에, 공간 필터 방법은 심장 모형을 재구성할 수 있다. 소스 벡터는 위치 r=(x,y,z) 및 시간 t 에서 이산 체적 공간(discrete volume space) V에서 s(r,t)로 정의된다. 소스들의 총 개수는 Q 이다. 상기 소스 벡터는 소스 전류의 크기와 방향을 나타낼 수 있다.
하나의 전류원으로부터 세 방향으로 계산된 리드필드벡터(l(r))는 수학식 2와 같이 나타낼 수 있고, 전체 전류원 Q에 대한 총 리드필드벡터(L(r))는 수학식 3과 같이 나타낼 수 있다(S130).
Figure 112014027564777-pat00003
Figure 112014027564777-pat00004
위치 r=(x,y,z)에서 시간이 t 일 때, 소스벡터를 s(r,t)로 정의하면, 심방과 심실의 소스벡터는 각각 sa(r,t)와 sb(r,t)로 나타낼 수 있다. 상기 리드필드벡터(L(r))와 소스벡터(sa(r,t), sb(r,t))를 선형결합하면 센서 공간에서 계산된 심방과 심실의 자기장(ba(t), bv(t))을 수학식 4와 같이 계산할 수 있다.
Figure 112014027564777-pat00005
여기서, V는 소스 공간을 나타낸다.
공간필터는 신호전송 및 수신을 위한 센서어레이(sensor array)에 사용되는 기술이다. 상기 공간필터(WT a,WT v)는 가중치벡터(weight vector)이다. 따라서, 측정된 자기장(심자도 신호)에 선형적으로 상기 공간필터(WT a,WT v)를 곱해줌으로써 소스전류의 신호를 추정할 수 있다. 즉, 상기 심방과 심실의 추정 소스전류 벡터는 수학식 5와 같이 나타낼 수 있다.
Figure 112014027564777-pat00006
상기 심방과 심실의 추정 소스전류 벡터(sa(r,t),sb(r,t))는 전체 추정 소스전류 벡터 s(r,t)에 포함되기 때문에, 추정 소스전류 벡터의 파워를 구하기 위해 일반적으로 유도되는 수학식에서는 심방 첨자(a)와 심실 첨자(v)로 표기하지 않고 설명한다. 상기 수학식 4와 수학식 5를 결합하면 추정 소스전류 벡터는 수학식 6과 같이 표현될 수 있다.
Figure 112014027564777-pat00007
여기서, r은 목표소스 위치이고, r'는 목표소스 이외의 위치이다.
Figure 112014027564777-pat00008
는 빔 응답(beam response)이고, 공간필터의 목표소스 위치가 아닌 다른 소스에서 발생된 왜곡된 누설전류(misleading leakage current)의 이득(gain)을 나타낸다.
이상적인 공간필터를 설계하기 위해서는 상기 빔 응답의 이득이 최소화되어야 한다. 따라서 최소크기기반(minimum-norm-based) 공간필터는 아래식의 제약을 조건으로 최적화 문제를 풀어 유도할 수 있다.
Figure 112014027564777-pat00009
Figure 112014027564777-pat00010
은 센서 어레이의 이득을 나타낸다. 수학식 7에서 C(W)는 비용함수(cost function)이고 전체 누설전류의 총합을 나타낸다. 상기 비용함수 C(W)는 수학식 8과 같이 나타낼 수 있다.
Figure 112014027564777-pat00011
tr{}은 행렬의 대각 성분의 합(trace)이고, δ(r)은 델타 함수, H(r')는 임의의 가중치행렬이다.
심방과 심실의 그램행렬(Gram matrix, Ga, Gv)은 상기 리드필드벡터와 심방과 심실의 임의의 가중치행렬을 결합하여 수학식 9와 같이 각각 정의할 수 있다.
Figure 112014027564777-pat00012
같은 크기를 가지는 소스전류에 대해, SQUID 센서에서 등가전류쌍극자의 거리가 가까운 경우 리드필드벡터의 크기는 크다. 반면 센서에서 먼 경우 리드필드벡터의 크기는 작다. 따라서 각각의 위치에 대한 리드필드벡터를 규격화(normalization)하면 멀리 있는 소스의 크기를 상대적으로 증가시킬 수 있다. 규격화된 리드필드벡터(
Figure 112014027564777-pat00013
)는 수학식 10과 같이 나타낼 수 있다(S140).
Figure 112014027564777-pat00014
상기 규격화된 리드필드벡터(
Figure 112014027564777-pat00015
)와 상기 그램행렬(G(r))을 선형결합하면 공간필터 또는 가중 벡터를 아래와 같이 구할 수 있다.
Figure 112014027564777-pat00016
최소분산(minimum variance;MV) 공간필터는 생체전자기(bioelectromagnetism) 분야에서 가장 많이 사용되는 적응공간필터이다. 공분산행은 측정데이터(b(t))의 앙상블 평균(ensemble average)으로 다음과 같이
Figure 112014027564777-pat00017
로 정의할 수 있다.
상기 최소분산 공간필터는 아래식의 제약을 조건으로 최적화 문제를 풀어 유도할 수 있다.
Figure 112014027564777-pat00018
최소분산 공간필터는 상기 규격화된 리드필드벡터(
Figure 112014027564777-pat00019
)와 상기 공분산행렬(D)의 선형 결합으로 아래와 같이 구할 수 있다.
Figure 112014027564777-pat00020
수학식 11과 13의 유일한 차이점은 그램행렬과 공분산행렬이다. 상기 공분산행렬(D)에 수학식 4를 적용하면 아래와 같이 나타낼 수 있다.
Figure 112014027564777-pat00021
수학식 9와 14를 비교해 보면, 임의의 가중치행렬 (H(r))이 소스 파워행렬
Figure 112014027564777-pat00022
과 유사함을 보여준다. 따라서 상기 가중치행렬(H(r))을 소스 파워행렬로 대체하면, 최소크기 공간필터의 성능이 최소분산 공간필터와 유사함을 기대할 수 있다.
더 정확한 결과를 얻기 위해서 그램 행렬에 노이즈 영향을 고려해야 한다. 함수의 조건수(condition number)는 입력변수의 작은 변화에 따른 출력값이 얼마나 변하는가를 나타내는 수이다. 실제 활용하는데 있어 그램 행렬의 조건수는 아주 크다. 구체적으로, 그램 행렬에서 노이즈 성분은 아주 작은 값이지만 그램 역행열로 변환되면 작은 노이즈값이 크게 변한다. 따라서, 그램 행렬의 역행렬을 계산하면 부정확한 결과가 초래된다. 정규화 방법(regularization method)을 적용하면 노이즈 성분을 제거함으로써 그램 역행렬 계산 에러를 줄일 수 있다. 따라서 심방과 심실의 수정된 배열이득제한최소크기 공간필터는 수학식 15와 같이 각각 표현할 수 있다.
Figure 112014027564777-pat00023
여기서 λa와 λv는 심방과 심실의 정규화 매개변수(regularization parameter)이고, 상기 정규화 매개변수는 측정값의 신호대 잡음비로 정해진다.
수학식 14에서 심방과 심실의 소스전류 벡터 sa(r,t)와 sv(r,t)는 알 수 없기 때문에 가중치행렬 Ha(r)과 Hv(r)을 구하기 위해 상기 소스전류 벡터는 추정 소스전류 벡터
Figure 112014027564777-pat00024
Figure 112014027564777-pat00025
로 대체할 수 있다(S193).
최적의 가중치행렬(H(r))을 구하기 위해 아래의 과정이 필요하다. 첫째, 가중치행렬을 단위행렬로 초기화한다(S150). 둘째, 수학식 15로부터 배열이득제한최소크기 공간필터의 출력값을 구한다(S180). 셋째, 수학식 5로부터 추정된 소스전류벡터(
Figure 112014027564777-pat00026
,
Figure 112014027564777-pat00027
)를 구한다(S191). 넷째 수학식 9의 가중치행렬(H(r))은 추정된 소스전류 벡터로 업데이트된다. 상기 가중치행렬(H(r))은 업데이트된 가중치가 일정한 값에 수렴할 때 까지 상기 첫째부터 넷째까지의 과정을 통해 회귀적으로 가중치행렬(H(r))은 업데이트된다(S193). 최종적으로 상기 업데이트가 완료되면 수학식 5로부터 추정 소스전류벡터가 구해질 수 있다.
상기 추정 소스전류 벡터(
Figure 112014027564777-pat00028
,
Figure 112014027564777-pat00029
)로부터 각각의 소스의 파워, 즉
Figure 112014027564777-pat00030
Figure 112014027564777-pat00031
를 구할 수 있다. 구해진 추정 소스 전류 파워는 최대 크기로 나누어 줌으로써 규격화(normalization)할 수 있다(S194).
상기 규격화된 추정 소스전류의 파워를 신호대잡음비로부터 설정한 제1 임계값(threshold)을 적용하여 심방과 심실의 제1 윤곽을 각각 재구성한다. 예를 들어, 제1 임계값이 0.5로 설정된 경우, 상기 제1 임계값 이상의 값을 가지는 복셀들만이 소스 공간에 표시될 수 있다. 제1 임계값 이상의 복셀들의 최외곽 복셀은 서로 연결되어 제1 윤곽을 제공할 수 있다. 제1-1 임계값은 심방에 대하여 설정되고, 제1-2 임계값은 심실에 대하여 설정될 수 있다.
제1 윤곽은 심방 및 심실에 대하여 각각 형성될 수 있다. 상기 심방과 심실의 제1 임계값들은 서로 다를 수 있다.
[결맞음(coherence) 매핑]
심근전류소스의 파워를 계산할 때 생리적 노이즈, 예를 들어 호흡이나 소화작용 등은 소스전류의 파워를 계산하는데 잘못된 영향을 미칠 수 있다. 따라서 더 정확한 심장모형을 생성하기 위해서 결맞음 매핑(coherence mapping) 방법을 추가로 적용할 수 있다. 결맞음은 두 개 이상의 파동이 위상에 따라 서로 간섭을 하는 현상을 의미한다. 심장의 경우, 영역에 따라 심근세포막의 이온투과율이 다르다. 따라서 근접 영역의 경우 심근 전류는 유사한 활동전위를 나타낸다. 하지만 먼 영역의 경우 상이한 활동전위(action potential)를 보여준다. 따라서 심방과 심실의 최대 신호원과 결맞는 신호원의 위치를 추정하면 심방과 심실의 영역으로 각각 윤곽을 재구성 할 수 있다.
도 6은 본 발명의 일 실시예에 따른 심방의 결맞음 매핑을 나타내는 도면이다.
도 7은 본 발명의 일 실시예에 따른 심실의 결맞음 매핑을 나타내는 도면이다.
도 6 및 도 7을 참조하면, 결맞음 분석은 주파수에 따른 신호간의 상관성을 나타내는 대표적인 분석법으로써 정상상태(stationary state) 두 신호에 대한 결맞음 함수는 상호 스펙트럼 밀도(cross spectral density)함수와 자기 스펙트럼 밀도(auto spectral density)함수로부터 계산할 수 있다.
소스 영역에서 추정된 신호파형(추정 소스전류 벡터)은 수학식 5로부터 공간필터와 측정된 자기장의 곱으로부터 재구성된다. Q번째 위치에서 재구성된 파형(추정 소스전류 벡터)의 시계열(time series)은
Figure 112014027564777-pat00032
로 나타난다. 심방과 심실의 소스전류의 파워가 최대인 지점에서의 파형(추정 소스벡터)은 각각
Figure 112014027564777-pat00033
Figure 112014027564777-pat00034
로 나타낼 수 있다(S210). 재구성된 시계열 파형(추정 소스벡터)
Figure 112014027564777-pat00035
,
Figure 112014027564777-pat00036
,
Figure 112014027564777-pat00037
는 고속 퓨리에 변환(fast Fourier transforms)에 의하여 Sa(f), Sv(f), Sq(f)로 변환될 수 있다(S230).
심방과 심실의 상호-스펙트럼(cross-spectrum) 및 Q번째 소스의 자기-스펙트럼(auto-spectrum)은 아래와 같이 정의할 수 있다(S240).
Figure 112014027564777-pat00038
여기서 < >는 평균을 의미하고, 윗첨자 *는 공액 복소수(complex conjugate)를 나타낸다.
심방과 심실의 결맞음은 상기 상호-스펙트럼과 자기-스펙트럼의 조합으로 수학식 17과 같이 나타낼 수 있다(S250).
Figure 112014027564777-pat00039
구체적으로, 심방의 결맞음(Cohaq(f))은 심방과 모든 소스 공간 사이에서 구한 상호-스펙트럼(Γaq(f))을 각각의 자기-스펙트럼(Γaa(f),Γqq(f) )의 곱의 제곱근으로 나누어 줌으로써 규격화할 수 있다. 마찬가지로 심실의 결맞음(Cohvq(f))도 상기 심방과 같이 규격화 할 수 있다. 심실의 결맞음(Cohvq(f))은 심실과 모든 소스 공간 사이에서 구한 상호-스펙트럼(Γvq(f))을 각각의 자기-스펙트럼(Γvv(f),Γqq(f))의 곱의 제곱근으로 나누어 줌으로써 규격화할 수 있다.
상기 규격화된 심방 결맞음 및 심실 결맞음은 신호대 잡음비로부터 정해진 각각의 제2 임계값을 기준으로 제2 윤곽을 재구성할 수 있다(S260). 예를 들어, 제2 임계값이 0.5로 설정된 경우, 상기 제2 임계값 이상의 값을 가지는 복셀들만이 소스 공간에 표시될 수 있다. 제2 임계값 이상의 복셀들의 최외곽 복셀은 서로 연결되어 제2 윤곽을 제공할 수 있다. 제2-1 임계값은 심방에 대하여 설정되고, 제2-2 임계값은 심실에 대하여 설정될 수 있다.
제2 윤곽은 심방 및 심실에 대하여 각각 형성될 수 있다. 상기 심방과 심실의 제2 임계값들은 서로 다를 수 있다.
심방과 심실의 제1 윤곽과 제2 윤곽을 결합하여 제3 윤곽을 재구성 할 수 있다(S300). 상기 재구성 방법은 심방과 심실의 상기 규격화된 소스전류의 파워와 상기 규격화된 결맞음을 곱해준 후 각각의 제3 임계값을 기준으로 제3 윤곽을 재구성할 수 있다. 제3 윤곽은 심방 및 심실에 대하여 각각 형성될 수 있다. 상기 심방과 심실의 제3 임계값들은 서로 다를 수 있다.
이하, 본 발명의 유효성을 증명하기 위한 수치 시뮬레이션(numerical simulation)과 모형 실험(phantom experiment) 및 수치 시뮬레이션과 모형 실험의 결과가 설명된다.
도 8은 본 발명의 일 실시예에 따른 수치 시뮬레이션을 위한 가상 소스 공간 및 가상 심장을 나타내는 도면이다.
도 8을 참조하면, 가상 심장(204,206)은 가상 소스 공간(202) 내에 위치한다. 상기 가상 심장은 가상 심방(204)과 가상 심실(206)을 포함한다.
상기 가상 소스 공간(202)은 xyz 직각 좌표계를 기준으로 직육면체일 수 있다. 상기 직각 좌표계에서 임의의 명치(xiphoid)위치가 원점(fiducial point)으로 정해진다. 상기 가상 소스 공간(202)은 10 mm 단위로 분할된 복셀(voxels)로 구성될 수 있다. 상기 가상 소스 공간(202)은 x축 방향으로 150 mm이고, y축 방향으로 170 mm이고, z축 방향으로 150 mm일 수 있다. 이에 따라, 총 복셀의 개수(Q)는 3314개 일 수 있다.
상기 가상 심방(204)에는 199개의 복셀이 할당되었고, 상기 가상 심실(206)에는 248개의 복셀이 할당되었다. 상기 가상 심장(204,206)에는 등가전류쌍극자(ECD)가 소스모델로 설정되었고, 신호의 형태는 사인파(sine wave)이다. 상기 심방의 등가전류쌍극자의 크기는 0.25 uAm이고, 심실의 등가전류쌍극자의 크기는 1 uAm이다. 그리고 신호대 잡음비 20 dB의 가우시안 백색 잡음(Gaussian white noise)이 추가되었다.
상기 심방 사인파의 중심주파수는 10 Hz이고, 심실 사인파의 중심주파수는 8 Hz이다. 또한 평면도체모델이 심장의 도체모델로 설정되었다. 심장의 활동 전위가 전파되는 것을 모사하기 위해 위치에 따라 사인파의 다른 위상이 할당되었다. 심장의 수축(systole)과 이완(diastole)을 나타내기 위해 심저(cardiac base)부터 심첨(cardiac apex), 즉 y축 방향으로 위상이 30 도씩 순차적으로 증가되었다. 또한 심장의 가로면(transverse plane), 즉 z축 방향으로 위상이 26 도씩 순차적으로 증가되었다.
도 9는 본 발명의 일 실시예에 따른 수치 시뮬레이션의 결과를 나타내는 도면이다.
도 9를 참조하면, 설정된 소스모델로부터 생성된 자기장을 센서 공간에서 측정하였다(S110). 상기 회귀생신그램행렬 배열이득제한최소크기 공간필터를 적용하여 심방과 심실 소스전류의 파워를 구하였다(S191). 상기 소스전류의 파워를 신호대 잡음비로 설정된 제1-1 임계값과 제1-2 임계값으로 필터링하여 심방과 심실의 제1 윤곽을 재구성하였다(S195). 상기 재구성된 신호파형으로부터 모든 소스영역에 대한 심방과 심실의 결맞음을 구하였다(S250). 상기 결맞음을 신호대 잡음비로 설정된 제2-1 임계값과 제2-2 임계값으로 필터링하여 심방과 심실의 제2 윤곽을 재구성하였다(S260). 마지막으로 상기 제1 윤곽과 제2 윤곽을 결합하여 심장의 제3 윤곽을 재구성하였다(S300). 그 결과, 상기 재구성된 심장의 윤곽이 가상의 심장 윤곽과 거의 일치함을 나타낸다.
도 10은 본 발명의 일 실시예에 따른 모형 실험을 위한 모형 심장(cardiac phantom)을 나타내는 도면이다.
도 10을 참조하면, 상기 모형심장(a)은 거위알(goose egg) 형태로 제작되었다. 길이 10 mm의 12개 전류 쌍극자가 상기 모형심장의 표면 5 mm 상에 위치한다. 상기 모형심장은 유리섬유(fiberglass)로 제작된 몸통모형(torso phantom)내에 위치한다. 상기 모형심장은 관상면(coronal plane)으로 45도 기울어지고, 사상면(sagittal plane)으로 30도 기울어진다. 상기 몸통모형내에 전기전도도(conductivity)가 약 0.16 S/m인 0.9 % 생리식염수(saline solution)를 채운다.
상기 모형 심장은 가상 소스 공간 내에 위치한다. 상기 가상 소스 공간은 상기 수치 시뮬레이션에 사용되었던 가상 소스 공간과 동일하다. 상기 몸통모형(b)의 표면에서 임의의 명치(xiphoid)위치가 원점(fiducial point)으로 정해진다. 상기 모형 심장에는 등가전류쌍극자(ECD)가 소스모델로 설정되었다. 상기 모형 심장 등가전류쌍극자의 크기는 100 uAm이고, 신호의 형태는 5 Hz와 10 Hz 사인파(sine wave) 조합이다. 또한 평면도체모델이 심장의 도체모델로 설정되었다.
도 11은 도 10의 모형 심장을 이용한 수치 시뮬레이션 결과를 나타내는 도면이다.
도 11을 참조하면, 상기 수치 시뮬레이션과 마찬가지로 소스전류의 파워와 결맞음 매핑 방법으로 심장의 제1 윤곽과 제2 윤곽을 구하였다. 상기 제1 윤곽과 제2 윤곽을 결합하여 심장의 제3 윤곽을 재구성하였다. 그 결과, 상기 재구성된 심장의 윤곽이 모형 심장 윤곽과 거의 일치함을 나타낸다.
도 12는 신호대 잡음비에 따른 정규화 매개변수의 비율을 나타내는 도면이다.
도 12를 참조하면, 정규화 매개변수의 비율은 상기 그램행렬의 최대특이값(maximum singular value)에 대한 비율이다. 신호대 잡음비가 클수록 정규화 매개변수의 비율이 작아진다. 또한 심방의 신호가 심실의 신호보다 작기 때문에 정규화 매개변수의 비율이 작게 나타낸다.
본 발명은 도면에 도시된 실시예를 참고로 설명되었으나 이는 예시적인 것에 불과하며, 본 기술 분야의 통상의 지식을 가진 자라면 이로부터 다양한 변형 및 균등한 다른 실시예가 가능하다는 점을 이해할 것이다. 따라서 본 발명의 진정한 기술적 보호 범위는 첨부된 특허청구범위의 기술적 사상에 의하여 정해져야 할 것이다.
101: 자기차폐실
122: 처리부
112: 심자도 센서

Claims (10)

  1. 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 제1 윤곽을 구하는 단계;
    결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계; 및
    상기 제1 윤곽과 상기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성하는 단계를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  2. 제1 항에 있어서,
    심장의 제1 윤곽을 구하는 단계는:
    심자도 장치를 이용하여 심근에서 발생된 자기 신호를 측정하는 단계;
    다채널로 측정된 심자도 신호에서 심방과 심실을 나타내는 P파와 T파 구간의 심장 자기 신호를 추출하는 단계;
    심장의 단위 소스전류 및 심장을 근사화한 평면체적도체모델 (horizontally layered volume conductor model)과 심자도 센서의 위치정보를 이용하여 자기장(리드필드 벡터)를 계산하는 단계; 및
    거리에 따른 신호원의 크기를 보상하기 위하여 상기 리드필드 벡터를 규격화(normalization)하는 단계를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  3. 제2 항에 있어서,
    심장의 제1 윤곽을 구하는 단계는:
    임의의 가중치행렬을 단위 행렬로 초기화하는 단계;
    상기 리드필드 벡터와 임의의 가중치행렬을 조합하여 그램행렬을 구하는 단계;
    상기 구해진 그램행렬에 측정신호의 신호대 잡음비를 고려하여 정규화 매개변수(regularization parameter)를 설정하는 단계;
    상기 정규화된 그램행렬과 규격화된 리드필드 벡터를 조합하여 배열이득제한최소크기 공간필터(array-gain constraint minimum-norm spatial filter)를 구하는 단계;
    상기 구해진 배열이득제한최소크기 공간필터와 측정된 심방과 심실 신호를 이용하여 추정 소스전류 벡터의 파워를 구하는 단계;
    상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하는지 확인하는 단계;
    상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하지 않는 경우 가중치행렬을 추정 소스벡터를 이용하여 갱신하는 단계;
    상기 구해진 추정 소스 전류 벡터의 파워를 규격화하는 단계; 및
    규격화된 소스전류벡터의 파워를 신호대잡음비로 설정된 제1 임계값을 적용하여 심방과 심실의 제1 윤곽을 재구성하는 단계들 중에서 적어도 하나를 더 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  4. 제1 항에 있어서,
    결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하는 단계는:
    심방과 심실에서 소스전류벡터의 파워가 가장 큰 지점을 각각 구하는 단계;
    측정된 심자도 신호와 회귀갱신그램행렬 배열이득제한최소크기 공간필터를 선형결합하여 모든 소스영역에서 신호파형(추정 소스전류 벡터)을 재구성하는 단계;
    상기 재구성된 신호파형(추정 소스전류 벡터)을 고속 퓨리에 변환(fast Fourier transform)하는 단계;
    상기 고속 퓨리에 변환된 신호를 이용하여 심방과 심실에서 모든 소스영역에 대한 상호-스펙트럼(cross-spectrum)과 자기-스펙트럼(auto-spectrum)을 구하는 단계;
    상기 구해진 상호-스펙트럼과 자기-스펙트럼을 조합하여 심방의 결맞음과 심실의 결맞음을 구하는 단계;
    상기 구해진 심방과 심실의 결맞음을 규격화하는 단계; 및
    상기 규격화된 결맞음을 신호대 잡음비로 설정된 제2 임계값을 적용하여 심방과 심실의 제2 윤곽을 재구성하는 단계들 중에서 적어도 하나를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  5. 회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 윤곽을 구하는 단계를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  6. 제5 항에 있어서,
    심장의 윤곽을 구하는 단계는:
    다채널로 측정된 심자도 신호에서 심방과 심실을 나타내는 P파와 T파 구간의 심장 자기 신호를 추출하는 단계;
    심장의 단위 소스전류 및 심장을 근사화한 평면체적도체모델 (horizontally layered volume conductor model)과 심자도 센서의 위치정보를 이용하여 자기장(리드필드 벡터)를 계산하는 단계; 및
    거리에 따른 신호원의 크기를 보상하기 위하여 상기 리드필드 벡터를 규격화(normalization)하는 단계를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  7. 제6 항에 있어서,
    임의의 가중치행렬을 단위 행렬로 초기화하는 단계;
    상기 리드필드 벡터와 임의의 가중치행렬을 조합하여 그램행렬을 구하는 단계;
    상기 구해진 그램행렬에 측정신호의 신호대 잡음비를 고려하여 정규화 매개변수(regularization parameter)를 추가하는 단계;
    상기 정규화된 그램행렬과 규격화된 리드필드 벡터를 조합하여 배열이득제한최소크기 공간필터(array-gain constraint minimum-norm spatial filter)를 구하는 단계;
    상기 구해진 배열이득제한최소크기 공간필터와 측정된 심방과 심실 신호를 이용하여 추정 소스전류 벡터의 파워를 구하는 단계;
    상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하는지 확인하는 단계;
    상기 구해진 소스전류 벡터의 파워가 일정한 값으로 수렴하지 않는 경우 가중치행렬을 추정 소스벡터를 이용하여 갱신하는 단계;
    상기 구해진 추정 소스 전류 벡터의 파워를 규격화하는 단계;
    규격화된 소스전류벡터의 파워를 신호대잡음비로 설정된 임계값을 적용하여 심방과 심실의 윤곽을 재구성하는 단계를 더 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  8. 심방과 심실에서 추정 소스전류벡터의 파워가 가장 큰 지점을 각각 구하는 단계;
    측정된 심자도 신호와 회귀갱신그램행렬 배열이득제한최소크기 공간필터를 선형결합하여 모든 소스영역에서 신호파형(추정 소스전류 벡터)을 재구성하는 단계;
    상기 재구성된 신호파형(추정 소스전류 벡터)을 고속 퓨리에 변환(fast Fourier transform)하는 단계;
    상기 고속 퓨리에 변환된 신호를 이용하여 심방과 심실에서 모든 소스영역에 대한 상호-스펙트럼(cross-spectrum)과 자기-스펙트럼(auto-spectrum)을 구하는 단계;
    상기 구해진 상호-스펙트럼과 자기-스펙트럼을 조합하여 심방의 결맞음과 심실의 결맞음을 구하는 단계;
    상기 구해진 심방과 심실의 결맞음을 규격화하는 단계; 및
    상기 규격화된 결맞음을 신호대 잡음비로 설정된 임계값을 적용하여 심방과 심실의 윤곽을 재구성하는 단계를 포함하는 것을 특징으로 하는 심장 윤곽 재구성 방법.
  9. 제1 항 내지 제8 항 중에서 어느 한 항의 방법을 컴퓨터에서 실행시키기 위한 프로그램을 기록한 기록매체.
  10. 심자도 센서 및 자기차폐실을 포함하고 심자도 신호를 측정하는 심자도 측정 장치; 및
    상기 심자도 신호를 처리하여 심장의 윤곽을 재구성하는 처리부를 포함하고,
    상기 처리부는:
    회귀갱신그램행렬 배열이득제한최소크기(array-gain constraint minimum-norm with recursively updated gram matrix; AGMN-RUG) 공간필터를 이용하여 구한 소스전류벡터의 파워(source current power)로부터 심장의 제1 윤곽을 구하고,
    결맞음 매핑(coherence mapping) 방법을 적용하여 심장의 제2 윤곽을 구하고, 그리고
    상기 제1 윤곽과 상기 제2 윤곽을 상호 결합하여 심장의 제3 윤곽을 구성하는 것을 특징으로 하는 심자도 측정 시스템.
KR1020140033636A 2014-03-21 2014-03-21 3차원 심장 윤곽 재구성 방법 KR101525485B1 (ko)

Priority Applications (3)

Application Number Priority Date Filing Date Title
KR1020140033636A KR101525485B1 (ko) 2014-03-21 2014-03-21 3차원 심장 윤곽 재구성 방법
PCT/KR2015/002569 WO2015142029A1 (ko) 2014-03-21 2015-03-17 3차원 심장 윤곽 재구성 방법
CN201580015416.3A CN106132288B (zh) 2014-03-21 2015-03-17 三维心脏轮廓重构方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
KR1020140033636A KR101525485B1 (ko) 2014-03-21 2014-03-21 3차원 심장 윤곽 재구성 방법

Publications (1)

Publication Number Publication Date
KR101525485B1 true KR101525485B1 (ko) 2015-06-03

Family

ID=53505268

Family Applications (1)

Application Number Title Priority Date Filing Date
KR1020140033636A KR101525485B1 (ko) 2014-03-21 2014-03-21 3차원 심장 윤곽 재구성 방법

Country Status (3)

Country Link
KR (1) KR101525485B1 (ko)
CN (1) CN106132288B (ko)
WO (1) WO2015142029A1 (ko)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116965822A (zh) * 2023-09-25 2023-10-31 合肥工业大学 心磁彩色空间圆图生成及波段时间识别方法及存储介质

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108470209B (zh) * 2018-03-27 2021-06-04 北京工业大学 一种基于格拉姆矩阵正则化的卷积神经网可视化方法
CN109579693A (zh) * 2018-11-26 2019-04-05 中国科学院上海光学精密机械研究所 一种基于互相干度最优的成像处理方法
CN115985505B (zh) * 2023-01-19 2023-12-12 北京未磁科技有限公司 多维度融合的心肌缺血辅助诊断模型及其构建方法

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20130045613A (ko) * 2011-10-26 2013-05-06 한국표준과학연구원 비침습적 심근 전기활동 매핑 방법

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR100562099B1 (ko) * 2003-08-28 2006-03-17 한국전기연구원 다채널 심자도 신호를 이용한 심자도의 특징점 및 특징구간의 검출방법
US8688192B2 (en) * 2011-01-31 2014-04-01 Seiko Epson Corporation High-resolution magnetocardiogram restoration for cardiac electric current localization
US9089274B2 (en) * 2011-01-31 2015-07-28 Seiko Epson Corporation Denoise MCG measurements

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20130045613A (ko) * 2011-10-26 2013-05-06 한국표준과학연구원 비침습적 심근 전기활동 매핑 방법

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116965822A (zh) * 2023-09-25 2023-10-31 合肥工业大学 心磁彩色空间圆图生成及波段时间识别方法及存储介质
CN116965822B (zh) * 2023-09-25 2023-12-29 合肥工业大学 心磁彩色空间圆图生成及波段时间识别方法及存储介质

Also Published As

Publication number Publication date
CN106132288B (zh) 2019-04-26
CN106132288A (zh) 2016-11-16
WO2015142029A1 (ko) 2015-09-24

Similar Documents

Publication Publication Date Title
US10201288B2 (en) Multi-electrode mapping system
US11013444B2 (en) Method and device for determining and presenting surface charge and dipole densities on cardiac walls
US8553956B2 (en) 3D current reconstruction from 2D dense MCG images
US9125581B2 (en) Continuous modeling for dipole localization from 2D MCG images with unknown depth
KR101310747B1 (ko) 비침습적 심근 전기활동 매핑 방법
CN108324263B (zh) 一种基于低秩稀疏约束的无创心脏电生理反演方法
US6856830B2 (en) Method and apparatus of three dimension electrocardiographic imaging
Pullan et al. The inverse problem of electrocardiography
Chinchapatnam et al. Model-based imaging of cardiac apparent conductivity and local conduction velocity for diagnosis and planning of therapy
US4961428A (en) Non-invasive method and apparatus for describing the electrical activity of the surface of an interior organ
US9089274B2 (en) Denoise MCG measurements
KR101525485B1 (ko) 3차원 심장 윤곽 재구성 방법
CN110811596B (zh) 基于低秩与稀疏约束和非局部全变分的无创心脏电位重建方法
CN110393522B (zh) 一种基于图总变分约束的无创心脏电生理反演方法
Fang et al. Noninvasive imaging of epicardial and endocardial potentials with low rank and sparsity constraints
Ha et al. Three-dimensional reconstruction of a cardiac outline by magnetocardiography
Yadan et al. An expert review of the inverse problem in electrocardiographic imaging for the non-invasive identification of atrial fibrillation drivers
El Houari et al. Investigating transmembrane current source formulation for solving the ECG inverse problem
Bhat et al. A Study of Imaging the Cardiac Activation Sequences in Electrocardiology
Svehlikova et al. Two Approaches for Inverse PVC Localization from Clinical ECG Data Using Heart Surface Potentials
Kramm et al. Constructional Features of a Multielectrode Electrocardiology Screening System
Lunttila et al. Regularization in cardiac source imaging
Tao Advanced Source Reconstruction and Volume Conductor Modeling for Fetal Magnetocardiography
Ahrens The inverse problem of magnetocardiography: activation time imaging with the extended Kalman filter and the unscented Kalman filter
Hanna et al. Imaging of cardiac electrical sources using a novel spatio-temporal map-based regularization method

Legal Events

Date Code Title Description
E701 Decision to grant or registration of patent right
GRNT Written decision to grant
FPAY Annual fee payment

Payment date: 20180418

Year of fee payment: 4

FPAY Annual fee payment

Payment date: 20190417

Year of fee payment: 5