Wykrywanie zniekształceń architektonicznych we wcześniejszych mammografiach poprzez analizę zorientowanych wzorców

Podsumowanie

Demonstrujemy metody wykrywania zniekształceń architektonicznych w poprzednich mammografach. Zorientowane struktury są analizowane za pomocą filtrów Gabora i portretów fazowych w celu wykrycia miejsc promieniowania wzorców tkankowych. Każde miejsce jest scharakteryzowane i sklasyfikowane za pomocą miar reprezentujących wzory spikulacyjne. Metody te powinny pomóc w wykrywaniu raka piersi.

Streszczenie

Demonstrujemy metody wykrywania zniekształceń architektonicznych w poprzednich mammografiach przypadków raka interwałowego na podstawie analizy orientacji wzorów tkanki piersi w mammografii. Stawiamy hipotezę, że zniekształcenia architektoniczne modyfikują normalną orientację wzorów tkanki piersi na obrazach mammograficznych przed powstaniem mas lub guzów. W początkowych krokach naszych metod zorientowane struktury w danej mammografii są analizowane za pomocą filtrów Gabora i portretów fazowych w celu wykrycia węzłopodobnych miejsc promieniujących lub przecinających się wzorów tkankowych. Każde wykryte miejsce jest następnie charakteryzowane za pomocą wartości węzła, wymiaru fraktalnego i miary dyspersji kątowej specjalnie zaprojektowanej do reprezentowania pikantnych wzorców związanych ze zniekształceniami architektonicznymi.

Nasze metody zostały przetestowane z bazą danych zawierającą 106 wcześniejszych mammografii 56 przypadków raka interwałowego i 52 mammografów 13 normalnych przypadków, korzystając z funkcji opracowanych do charakterystyki zniekształceń architektonicznych, klasyfikacji wzorców za pomocą kwadratowej analizy dyskryminacyjnej i walidacji z procedurą pominięcia jednego pacjenta. Zgodnie z wynikami analizy charakterystyki operacyjnej odbiornika swobodnej odpowiedzi, nasze metody wykazały zdolność do wykrywania zniekształceń architektonicznych we wcześniejszych mammografiach, wykonanych średnio 15 miesięcy przed klinicznym rozpoznaniem raka piersi, z czułością 80% przy około pięciu fałszywie dodatnich wynikach na pacjentkę.

Wprowadzenie

Rak piersi jest główną chorobą dotykającą kobiety i jest drugą najczęstszą przyczyną zgonów związanych z rakiem wśród kobiet ^1,2. Aby zwiększyć szansę na przeżycie i rokowanie chorych pacjentek poprzez skuteczne leczenie we wczesnych stadiach raka piersi, choroba musi zostać wykryta tak wcześnie, jak to możliwe. W retrospektywnej analizie przypadków raka piersi zaobserwowano subtelne oznaki nieprawidłowości na wcześniej nabytych mammografiach przesiewowych ^3,4. Zniekształcenie architektoniczne jest jednym z takich zlokalizowanych objawów mammograficznych prawdopodobnie wczesnych stadiów raka piersi, który jest trudny do wykrycia ^5,6. Związane z tym wzorce są niejasno opisane jako zniekształcenie normalnej architektury piersi bez widocznej określonej masy. Zniekształcenia architektoniczne mogą pojawić się na początkowych etapach tworzenia się masy piersi lub guza. Stawiamy hipotezę, że mammografia przesiewowa uzyskana przed wykryciem raka piersi może zawierać subtelne oznaki wczesnych stadiów raka piersi, w szczególności zniekształcenia architektoniczne.

Rysunek 1a pokazuje wcześniejszy obraz mammograficzny przypadku raka wykrytego na ekranie. Obszar nieprawidłowości zidentyfikowany przez radiologa (J.E.L.D.) jest zaznaczony czerwonym prostokątem. Poprzednia mammografia została wykonana 24 miesiące przed mammografią detekcyjną pokazaną na rycinie 1b. Poprzednia mammografia została uznana za wolną od objawów raka podczas pierwotnego badania przesiewowego. W analizie retrospektywnej i w porównaniu z mammografią wykrywającą, podejrzany obszar związany z miejscem wykrycia raka został oznaczony przez radiologa i jest zaznaczony na czerwono na poprzedniej mammografii. W podejrzanym regionie znajdują się ślady zniekształceń architektonicznych, w tym kolce.

Techniki i systemy diagnostyki wspomaganej komputerowo (CAD) oferują potencjał do osiągnięcia zwiększonej czułości w wykrywaniu raka piersi ^2,7-9. Jednak w porównaniu z liczbą publikacji, które istnieją w literaturze na temat wykrywania innych objawów raka piersi, takich jak masy i zwapnienia, tylko niewielka liczba badań została zgłoszona na temat wykrywania zniekształceń architektonicznych przy braku centralnej masy ^10-17. Stwierdzono, że dostępne na rynku systemy CAD słabo radzą sobie z wykrywaniem zniekształceń architektonicznych ¹⁸. Badania nad wykrywaniem zniekształceń architektonicznych we wcześniejszych mammografiach u przypadków raka wykrytych w badaniach przesiewowych lub w odstępach międzynowotworowych ^3,4,19-22 mogą pomóc w opracowaniu strategii wykrywania i leczenia chorób piersi we wczesnych stadiach i prowadzić do poprawy rokowania dla pacjentki ²³.

Przygotowanie obrazów do eksperymentu

Eksperymenty przeprowadzono na 158 obrazach mammograficznych, w tym 106 wcześniejszych mammografiach 56 osób, u których zdiagnozowano raka piersi i 52 zdjęciach 13 normalnych osób. Zgodę na badanie uzyskano od Rady Etyki Badań nad Zdrowiem Conjoint, Biura Bioetyki Medycznej Uniwersytetu w Calgary oraz Regionalnego Urzędu ds. Zdrowia w Calgary. Obrazy uzyskano z Screen Test: Alberta Program for the Early Detection of Breast Cancer ^21,24,25.

Mammografia uzyskana podczas ostatniej zaplanowanej wizyty w programie badań przesiewowych przed diagnozą raka poza programem badań przesiewowych została oznaczona jako wcześniejsze mammografie przypadków raka interwałowego. Odpowiednie mammografie diagnostyczne nie były dostępne. Wszystkie ze 106 wcześniejszych mammografii, z wyjątkiem dwóch, zostały uznane za wolne od jakichkolwiek objawów raka piersi w momencie ich nabycia i analizy w programie badań przesiewowych; Osoby, które otrzymały pozostałe dwie mammografie, zostały skierowane na biopsję. Odstęp czasu między rozpoznaniem raka a wcześniejszą mammografią wynosił od 1,5 miesiąca do 24,5 miesiąca, przy średniej 15 miesięcy i odchyleniu standardowym wynoszącym 7 miesięcy. Wszystkie wcześniejsze mammografie przypadków raka interwałowego dostępne w bazie danych zostały włączone do niniejszego badania, z wyjątkiem sześciu obrazów, na których nie można było zidentyfikować żadnych podejrzanych części.

Mammografia ekranowo-filmowa została zdigitalizowana w rozdzielczości przestrzennej 50 μm i rozdzielczości w skali szarości 12 bitów na piksel przy użyciu skanera laserowego Lumiscan 85 (Lumisys, Sunnyvale, CA). Biegły radiolog specjalizujący się w mammografii (J.E.L.D.) dokonał przeglądu wszystkich 106 wcześniejszych mammografii przypadków raka interwałowego i oznaczył podejrzane obszary zniekształceń architektonicznych prostokątnymi pudełkami na podstawie dostępnych raportów z późniejszego obrazowania lub biopsji lub szczegółowej kontroli mammografii. Spośród 106 wcześniejszych obrazów mammograficznych w zbiorze danych wykorzystanym w niniejszym badaniu, 38 obrazów ma widoczne zniekształcenia architektoniczne, a pozostałe 68 obrazów zawiera wątpliwe lub nie zawiera wyraźnie widocznych zniekształceń architektonicznych. Każda wcześniejsza mammografia zawiera pojedyncze miejsce zniekształcenia architektonicznego, które zostało zidentyfikowane przez prostokątną ramkę narysowaną przez radiologa. Średnia szerokość, wysokość i obszar 106 podejrzanych części obrazów oznaczonych przez radiologa to 56 mm, 39 mm i 2274 mm², z odchyleniem standardowym wynoszącym odpowiednio 11,8 mm, 11,6 mm i 1073,9 mm².

Protokół

1. Przegląd metodologii

W naszej procedurze, potencjalne miejsca zniekształceń architektonicznych w mammografach są automatycznie wykrywane poprzez analizę zorientowanych wzorów teksturalnych z zastosowaniem banku filtrów Gabora ²⁶ i modelowanie portretów fazowych ^11,27. Wykryte miejsca są następnie przetwarzane przez etapy ekstrakcji cech lub miar w celu scharakteryzowania zniekształceń architektonicznych, opracowania wyszkolonego klasyfikatora i zastosowania algorytmu do rozpoznawania lub klasyfikacji wzorców. Procedura jest podsumowana następującymi krokami ^11,20,21:

1. Segmentuj część piersi na podanym obrazie mammograficznym za pomocą adaptacyjnego progowania i otworu morfologicznego.
2. Zastosuj zestaw 180 filtrów Gabora z kątami rozmieszczonymi równomiernie w zakresie od -90° do +90°, aby uzyskać obraz wielkości Gabora, M (i, j) i obraz kąta Gabora, θ (i, j), wybierając reakcję i kąt filtra o najwyższej odpowiedzi w każdym pikselu (i, j).
3. Wybierz interesujące nas struktury krzywoliniowe (CLS), takie jak kolce i tkanka włóknisto-gruczołowa, odróżniając je od struktur zakłócających, takich jak krawędzie mięśnia piersiowego, tkanka miąższowa, granica piersi i hałas, używając pola orientacji, pola gradientu, techniki niemaksymalnego tłumienia (NMS) i dodatkowych warunków ¹¹.
4. Przefiltrować pole orientacji za pomocą filtra Gaussa z odchyleniem standardowym wynoszącym 7 pikseli i czterokrotnie zmniejszyć próbkowanie w celu zmniejszenia szumu i dalszych wymagań obliczeniowych ^11,20.
5. Zastosuj liniowe modelowanie fazowo-pionowe, z przesuwanym oknem analizy o rozmiarze 10 x 10 pikseli przy 800 μm/piksel, z jednym pikselem na krok, do filtrowanego pola orientacji, z określonymi warunkami wyboru map portretu fazowego związanych z określonymi typami wzorców węzłów ^11,20.
6. Oddaj głos, jeśli spełnione są określone warunki, w pozycji wyznaczonej przez stały punkt dla każdej pozycji okna analizy, aby utworzyć mapę węzłów.
7. Przefiltruj mapę węzłów za pomocą okna Gaussa o rozmiarze 35 x 35 pikseli, z empirycznie określonym odchyleniem standardowym wynoszącym 6 pikseli (4,8 mm), aby skonsolidować głosy znajdujące się w bliskiej odległości od siebie.
8. Przeanalizuj mapę węzłów, porządkując szczyty na mapie węzłów.
9. Wytnij obszary zainteresowania (ROI) o rozmiarze 128 x 128 pikseli, z wyjątkiem krawędzi, od oryginalnego obrazu, przy czym środek każdego obszaru zainteresowania znajduje się w środku powiązanego piku na mapie węzłów. Na krawędziach przetwarzanego obrazu utwórz ROI, aby uwzględnić tyle danych obrazu, ile jest dostępnych w określonym oknie.
10. Wyprowadź cechy lub miary, aby scharakteryzować wzory spikulacyjne związane ze zniekształceniem architektonicznym i oddziel je od normalnych wzorców tkankowych, które spełniły niektóre z warunków początkowych.
11. Opracuj wyszkolony klasyfikator do rozróżniania cech miejsc ze zniekształceniami architektonicznymi od cech normalnych wzorców tkanek przy użyciu zestawu treningowego ROI sklasyfikowanych przez radiologa.
12. Zastosuj przeszkolony klasyfikator do zestawu przypadków testowych i zweryfikuj wyniki z diagnozą dostarczoną przez radiologa i opartą na biopsji.

Kroki 1-9 wymienione powyżej są automatycznie stosowane do danego obrazu mammograficznego. Wybrane kroki procedury wymienione powyżej zostały opisane i zilustrowane w poniższych sekcjach.

2. Wstępne przetwarzanie obrazów mammograficznych

Etap wstępnego przetwarzania składa się z następujących kroków:

1. Przefiltruj podany obraz mammograficzny za pomocą filtra Gaussa, o odchyleniu standardowym 2 pikseli i rozmiarze 13 x 13 pikseli przy rozdzielczości 50 μm/piksel i 12 bitach/piksel, a następnie zmniejsz próbkowanie do 200 μm/piksel i rozdzielczości 8 bitów/piksel
.
2. Odbij obraz, jeśli przedstawia prawą pierś.
3. Segmentuj obszar piersi na obrazie mammograficznym przy użyciu adaptacyjnej metody progowania Otsu i otwarcia morfologicznego za pomocą elementu strukturyzującego w kształcie dysku o promieniu 25 pikseli (5 mm przy 200 μm / piksel) ^21,28,29.
4. Wykryj przybliżoną granicę piersi ^10,21.

Rysunek 2A pokazuje oryginalną wcześniejszą mammografię. Rysunek 2B tego samego rysunku przedstawia wynik przybliżonej segmentacji części piersiowej, która jest wykorzystywana w kolejnych etapach przetwarzania i analizy.

3. Ekstrakcja zorientowanych wzorców za pomocą filtrów Gabora

Rzeczywista funkcja filtra Gabora zorientowana na -90° jest określona w naszej pracy jako ^10,30:

gdzie σ_x i σ_y są wartościami odchylenia standardowego w kierunkach x i y, a ƒ_o jest częstotliwością modulującej sinusoidy. Filtry pod innymi kątami uzyskuje się przez obrócenie tej funkcji za pomocą transformacji współrzędnych jako:

gdzie (x', y') jest zbiorem współrzędnych obróconych o kąt α.

Parametry w równaniu 1 do filtrowania mammografii są uzyskiwane w naszej pracy poprzez uwzględnienie średniej wielkości wzorów tkanki piersi, które mają być wykryte, jak następuje: ¹⁰:

1. Niech Τ będzie pełną szerokością w połowie maksimum składnika Gaussa w równaniu 1 wzdłuż osi x.
2. Niech Τ = 4 piksele, co odpowiada grubości 0,8 mm przy rozmiarze piksela 200 μm.
3. Oblicz .
4. Niech okres członu cosinusa będzie równy Τ, a następnie ƒ_o = 1/Τ.
5. Niech wartość σ_y zostanie zdefiniowana jako σ_y = lσ_x, gdzie l określa wydłużenie filtra Gabora w kierunku y, w porównaniu z szerokością filtra w kierunku x. Do analizy mammografii przy 200 μm/piksel należy użyć l = 8.

Bank 180 prawdziwych filtrów Gabora, równomiernie rozmieszczonych w zakresie od -90° do +90°, jest używany w naszych metodach do wykrywania zorientowanych wzorców w mammografach ^{10, 21}. Dla każdego podanego obrazu obraz Gabora (M(i, j) i obraz kąta Gabora (θ(i, j) są uzyskiwane przy użyciu odpowiedzi i kąta filtra Gabora z najwyższą odpowiedzią dla każdego piksela (i, j).

Filtr Gabora ma niezerową odpowiedź wielkości na początku płaszczyzny częstotliwości (częstotliwość zerowa). Ponieważ składowe o niskiej częstotliwości nie są związane z obecnością zniekształceń architektonicznych, pożądane jest zmniejszenie wpływu składowych o niskiej częstotliwości obrazu mammograficznego na wielkość pola orientacyjnego. W związku z tym obrazy mammograficzne są filtrowane górnoprzepustowo przed ekstrakcją pola orientacyjnego. Osiąga się to poprzez obliczenie różnicy między oryginalnym obrazem a wersją tego samego obrazu z filtrem dolnoprzepustowym. Filtr dolnoprzepustowy używany w tym kroku jest filtrem Gaussa o odchyleniu standardowym równym σ_y zdefiniowanym jak powyżej.

Chociaż można by zapisać przefiltrowany obraz dla każdego interesującego kąta, w obecnej pracy, maksymalna odpowiedź na każdy piksel w stosunku do wszystkich użytych filtrów (kątów) jest zapisywana w pojedynczym obrazie, określana jako odpowiedź wielkości Gabora; odpowiedni kąt filtra Gabora jest zapisywany w każdym pikselu na innym obrazie, określana jako odpowiedź kąta Gabora. Razem dwa obrazy wyjściowe zapewniają pole orientacji danego obrazu.

Rysunek 3A przedstawia testowy obraz rośliny. Rysunek 3B przedstawia widmo Fouriera obrazu, które przedstawia koncentrację energii pod różnymi kątami. Wszystkie części obrazu o tej samej orientacji, niezależnie od ich położenia i rozmiaru, mają składowe częstotliwości (energia widmowa) zlokalizowane w paśmie kątowym lub sektorze umieszczonym pod kątem 90° w stosunku do ich orientacji na obrazie. Wyniki filtrowania obrazu za pomocą filtrów Gabora o Τ = 8 pikseli i l = 8 przedstawiono na rysunkach 3C i D. Oczywiste jest, że filtry Gabora wyodrębniły części rośliny zorientowane pod różnymi kątami z dużą reakcją na wielkość i że reakcja kątowa zgadza się z orientacją dominującej cechy obecnej w odpowiednim pikselu. Korzystając z zestawu filtrów Gabora zorientowanych pod kilkoma kątami w zakresie od -90° do +90°, wyodrębniliśmy wszystkie zorientowane komponenty obecne na obrazie i ich kąty w każdym pikselu. Oczywiste jest, że reakcja filtrów Gabora jest prawie zerowa w gładkich obszarach o tym samym poziomie intensywności i bez struktur o preferowanej orientacji, takich jak części doniczki i ściana.

Rysunek 4 pokazuje odpowiedzi na wielkość Gabora i kąt uzyskane dla mammografii ze zniekształceniem architektonicznym pokazanym na rysunku 2B. Oczywiste jest, że filtry Gabor wyodrębniły zorientowane komponenty o wysokiej reakcji, a także odpowiednie kąty. Widać również, że reakcja filtrów Gabora jest niska w gładkich obszarach o prawie stałej gęstości i bez struktur o preferowanej orientacji. Po bliższym przyjrzeniu się można również zauważyć, że reakcja filtrów Gabora zależy od kontrastu zorientowanej struktury w stosunku do jej tła, a nie tylko od jej gęstości czy jasności. Wyniki te wynikają z pasmowoprzepustowego charakteru filtrów Gabora.

4. Wybór struktur krzywoliniowych

Mammografie zawierają wiele CLS odpowiadających przewodom, naczyniom, więzadłom, tkance miąższowej i krawędziom mięśnia piersiowego. Niektóre nieprawidłowości w mammografii mogą charakteryzować się obecnością pewnych typów CLS, takich jak spikulowane masy ^12,31,32 i zniekształcenia architektoniczne ^10,11,33, lub asymetryczną strukturą zorientowanej tekstury w obrazie piersi ³⁴. Z drugiej strony, niektóre rodzaje zmian, takie jak ograniczone masy, mogą być zasłonięte przez kilka CLS nałożonych na zmiany w wyświetlanych obrazach mammograficznych; Wygląd takich zmian może być zmieniony i może prowadzić do fałszywie ujemnego wykrycia lub błędnej diagnozy. Analiza CLS obecnych w mammografach może poprawić wydajność algorytmów wykrywania mas spikulowanych i zniekształceń architektonicznych, jak sugerują Zwiggelaar i wsp. ³⁵. W związku z tym identyfikacja CLS jest ważnym krokiem w wykrywaniu zniekształceń architektonicznych.

Chociaż bank filtrów Gabora używany w niniejszej pracy jest wrażliwy na struktury liniowe, takie jak kolce i włókna, wykrywa również inne silne krawędzie, takie jak krawędzie mięśnia piersiowego, krawędzie tkanki miąższowej i ściany naczyń, jako struktury zorientowane. Mocne krawędzie wokół dysku włóknisto-gruczołowego ³⁶ mogą być wykorzystane do wykrywania szczególnej formy zniekształcenia architektonicznego ³⁷ znanego jako retrakcja ogniskowa. Jednak w niniejszej pracy ważne jest, aby tylko CLS związane z tkankami włóknisto-gruczołowymi były identyfikowane jako cechy zorientowane.

Metoda wyboru CLS-ów w obecnej pracy obejmuje następujące trzy kroki:

1. Na danej mammografii należy podzielić obszar piersi na segmenty zgodnie z opisem w punkcie 2.
2. Wykrywanie rdzeniowych pikseli CLS poprzez zastosowanie techniki NMS ^35,38 do obrazu odpowiedzi Gabora na magnitudę.
3. Odrzuć piksele CLS w miejscach o silnym nachyleniu ³³.

Algorytm NMS identyfikuje kluczowe piksele CLS, porównując każdy piksel na obrazie odpowiedzi magnitudy z jego sąsiadami wzdłuż kierunku, który jest prostopadły do lokalnego kąta pola orientacji; patrz Rysunek 5. Jeśli badany piksel ma większą wartość wielkości niż odpowiadający mu sąsiadowie, piksel jest podstawowym pikselem CLS. NMS jest częstym krokiem w wielu detektorach krawędziowych (takich jak detektor krawędzi Canny ³⁹). 35 Zwiggelaar i in. ³⁵ wykorzystali NMS do wykrywania pikseli CLS w taki sam sposób, jak opisano w tej sekcji.

Obecność silnego gradientu może spowodować falowanie w odpowiedzi na wielkość Gabora, prowadząc do błędnego wykrycia CLS. Rdzeniowe piksele CLS związane z obecnością silnych gradientów są odrzucane przez kryteria zaproponowane przez Karssemeijera i te Brake ¹² w kontekście wykrywania zmian kolczastych. Gradient obrazu mammograficznego uzyskuje się przy użyciu pierwszej pochodnej Gaussa o odchyleniu standardowym wynoszącym pięć pikseli (1 mm). Dla każdego rdzenia piksela CLS kierunek gradientu jest porównywany z kierunkiem pola orientacji. Jeśli różnica między kierunkiem pola orientacji a kierunkiem prostopadłym do gradientu jest mniejsza niż 30°, odpowiadający mu piksel CLS rdzenia jest odrzucany.

CLS w krążku włóknisto-gruczołowym zazwyczaj mają zmniejszony kontrast w porównaniu do CLS na zewnątrz dysku włóknisto-gruczołowego. W związku z tym CLS w krążku włóknisto-gruczołowym mają mniejsze wartości odpowiedzi wielkości Gabora niż CLS na zewnątrz dysku. W celu przypisania takiej samej wagi wszystkim pikselom CLS niezależnie od lokalizacji oraz w celu zapewnienia wykrycia odpowiednich pikseli CLS o niskim kontraście, takich jak kolce w krążku włóknisto-gruczołowym, pole wielkości M (i, j) zastępuje się do dalszego przetwarzania obrazem złożonym wyłącznie z rdzeniowych pikseli _{CLS, M CLS(}i, j), zdefiniowane w następujący sposób:

Obraz M_CLS(i, j) przekazuje ważne informacje na temat obecności CLS. Rysunek 6 przedstawia wyniki selekcji CLS z pełną mammografią i ROI. Ponieważ na obecność zniekształceń architektonicznych wskazuje geometryczny układ powiązanych CLS, a nie ich gęstość lub intensywność, wielkość wykrytych CLS ma mniejsze znaczenie niż układ przestrzenny zorientowanych struktur.

5. Wykrywanie i oznaczanie podejrzanych stron poprzez analizę portretów fazowych

Rao i Jain ⁴⁰ opracowali metodę analizy zorientowanej tekstury na obrazach poprzez powiązanie odpowiedniego pola orientacji gradientu z wyglądem portretów fazowych. Portret fazowy układu dwóch równań różniczkowych pierwszego rzędu pokazuje możliwe trajektorie zmiennych stanu ²⁷.

Niech p(t) i q(t), t R, reprezentują dwie różniczkowalne funkcje czasu t, powiązane jako

Tutaj p^•(t) i q^•(t) są pochodnymi pierwszego rzędu względem czasu, a F i G są funkcjami p i q ¹⁰. Biorąc pod uwagę warunki początkowe p(0) i q(0), roztwór [p(t), q(t)] można przedstawić w postaci parametrycznej trajektorii lub strumienia hipotetycznej cząstki w płaszczyźnie (p, q). Cząstka znajduje się w [p(0), q(0)] w czasie t = 0 i porusza się po płaszczyźnie (p, q) z prędkością [p^•(t) i q^•(t)]. Płaszczyzna (p, q) jest znana jako płaszczyzna fazowa układu. Portret fazowy to wykres możliwych trajektorii cząstki w płaszczyźnie fazowej. Punkt stały to punkt na płaszczyźnie fazy, w którym p^•(t) = 0 i q^•(t) = 0. Cząstka pozostawiona w ustalonym punkcie pozostaje nieruchoma. W przypadku układu afinicznego mamy

Tutaj A jest macierzą 2 x 2, a b jest macierzą kolumnową 2 x 1. Środek (p₀ , q₀) portretu fazowego jest określony przez punkt stały jako

Jeśli skojarzymy funkcje p(t) i q(t) ze współrzędnymi x i y płaszczyzny przetwarzanego obrazu, odpowiednim polem orientacji jest

Tutaj Φ(x y) jest kątem wektora prędkości [ p^•(t), q^•(t)] względem osi x w (x, y) = [p(t), q(t)]. Kojarzymy Φ(x y) z reakcją kąta Gabora θ (i, j) i definiujemy funkcję błędu, która ma być zminimalizowana jako

gdzie [a, b] i [c, d] to dwa rzędy A. Ostatni warunek przewiduje wyższą karę (koszt) za odchylenie w konfiguracjach macierzy A od tych związanych z wzorcami węzłów spikulowanych. Równanie podane powyżej reprezentuje Φ(x, y) na dyskretnej siatce (i, j) zamiast w przestrzeni ciągłej (x, y). Oszacowania A i b, które minimalizują ε²(A, b), uzyskuje się za pomocą następującej procedury:

1. Uzyskaj wstępne oszacowania A i b poprzez minimalizację ε²(A, b) przy użyciu metody symulowanego wyżarzania ⁴¹.
2. Uzyskanie optymalnych oszacowań przez uściślinie początkowych oszacowań przy użyciu nieliniowego algorytmu najmniejszych kwadratów ⁴².

W modelu opisanym powyżej, istnieją trzy możliwe typy portretów fazowych: węzł, siodło i spirala. Typ portretu fazowego jest określany przez wartości własne A ^10,27,30,40. Pole orientacji obrazu z teksturą można opisać, określając typ portretu fazowego, który jest najbardziej podobny do jego pola orientacji. Ponieważ wzory spiralne nie są interesujące w analizie mammografii, ograniczamy matrycę A do symetrii, co skutkuje tylko dwoma typami portretów fazowych: węzłem i siodłem.

Ze względu na spodziewaną obecność wielu kolców pod różnymi kątami, które nakładają się na wyświetlany obraz mammograficzny, przypuszczamy, że miejsce zniekształceń architektonicznych będzie miało cechy podobne do węzłów. Jednak normalne tkanki, przewody, naczynia i inne zorientowane struktury w piersi mogą również zostać rzutowane i nałożone, tworząc wzory, które naśladują pojawienie się zniekształceń architektonicznych na mammografii. W związku z tym analizujemy mapę węzłów pod kątem wykrywania podejrzanych lokalizacji lub potencjalnych miejsc zniekształceń architektonicznych i analizujemy wykryte miejsca poprzez dalsze kroki ekstrakcji cech i klasyfikacji wzorców.

Ponieważ mammografia może wykazywać kilka wzorców, stosujemy przesuwane okno analizy o rozmiarze 10 x 10 pikseli, przy 800 μm/piksel, z jednym pikselem na krok. Dla każdej pozycji okna oddawany jest głos na mapie, zwanej mapą węzłów, w pozycji wyznaczonej przez odpowiadający jej punkt stały, jeśli spełnione są wszystkie zastosowane warunki. Wyniki związane z macierzą A o numerze warunku większym niż 3,0 są odrzucane w celu zignorowania wzorców, które nie powinny być związane ze zniekształceniem architektonicznym ¹¹. Ponadto na odległość między punktem stałym a położeniem odpowiedniego okna analizy nałożony jest dodatkowy warunek: jeśli odległość jest mniejsza niż trzy piksele (2,4 mm) lub większa niż 20 pikseli (16 mm), wyniki dla bieżącego okna analizy są odrzucane. Wielkość głosu ustala się jako równy stosunkowi miary dopasowania ε² (A, b), zdefiniowanej w równaniu 7, do liczby warunku A, aby podkreślić izotropię portretu fazowego. Mapa węzłów jest następnie analizowana w celu wykrycia lokalnych maksimów lub szczytów, które powinny wskazywać miejsca zniekształceń architektonicznych. Jednak procedura ta powoduje również wykrycie wielu miejsc fałszywie dodatnich (FP) z powodu nałożonych na siebie normalnych struktur.

Przy każdym szczycie na mapie węzłów, automatycznie wyodrębniamy ROI, o rozmiarze 128 x 128 pikseli, z wyjątkiem krawędzi obrazów, z obrazu mammograficznego przy 200 μm/piksel. Oznaczamy ROI w miejscach wskazanych przez piki na mapie węzłów, w kolejności malejącej wartości szczytów, aż do maksymalnej liczby 30 ROI na mammografię.

Gdy mammografia ze znanymi diagnozami jest używana do trenowania naszej procedury, automatycznie wykryte ROI wraz z ich centrami w częściach zniekształceń architektonicznych zidentyfikowanych przez radiologa są oznaczane jako prawdziwie dodatnie (TP) ROI; inne są oznaczane jako ROI FP do wykorzystania w procedurze treningowej. Gdy mammografia jest analizowana przy użyciu wytrenowanej procedury, wszystkie ROI wykryte jak powyżej są przetwarzane w celu klasyfikacji bez żadnego oznakowania.

Rysunek 7 pokazuje mapę węzłów i wykryte ROI dla mammografii pokazanej na rysunku 2B. Czerwone prostokąty wskazują podejrzany obszar zaznaczony przez radiologa.

Rysunek 8 pokazuje liczbę ROI TP i FP uzyskanych automatycznie z kilku mammografii. Większość TP ROI ma kilka spikul i zorientowanych wzorów rozmieszczonych w szerokim zakresie kątów. Z drugiej strony ROI FP mają mniejszą liczbę normalnych struktur tkankowych zorientowanych w węższym zakresie kątów; Niezależnie od tego, ze względu na ich superpozycję w rzutowanym obrazie mammograficznym, naśladują one węzłowe cechy zniekształceń architektonicznych.

Nasza strategia polega na wykrywaniu podejrzanych regionów z wysoką skutecznością lub czułością na początkowym etapie (z odpowiednio niską liczbą wyników fałszywie ujemnych), nawet jeśli towarzysząca im liczba FP jest duża. Kolejny etap analizy ROI ma na celu pomoc w zmniejszeniu FP poprzez efektywną charakterystykę i klasyfikację wykrytych ROI.

6. Charakterystyka zniekształceń architektonicznych

Automatycznie wykryty ROI, w tym zniekształcenia architektoniczne, wyśrodkowany na szczycie na mapie powiązanych węzłów, prawdopodobnie posiada kilka pikulców rozrzuconych pod różnymi kątami. Spodziewamy się, że ta cecha doprowadzi do szerokiego kątowego rozrzutu energii w domenie obrazu i energii widmowej w domenie Fouriera. W naszych poprzednich pracach wykazaliśmy, że taką dyspersję kątową można skutecznie przedstawić w postaci diagramu różanego, który jest histogramem kątowym ^21,22. Normalizujemy diagram róży tak, aby miał jednostkę powierzchni i traktujemy go jako funkcję gęstości prawdopodobieństwa (PDF). Następnie charakteryzujemy plik PDF każdego zwrotu z inwestycji za pomocą entropii, która jest statystyczną miarą nieporządku lub rozproszenia.

Zwiększone rozproszenie wzorów tkanek w regionach z zniekształceniami architektonicznymi modyfikuje fraktalną naturę normalnej tkanki piersi. Powszechnie stosowane modele fraktali opierają się na wieloskalowych zagnieżdżonych wzorcach samopodobnych ^43-46. Innym modelem zachowania fraktalnego jest ułamkowy ruch Browna (fBm), który jest związany z widmem w dziedzinie częstotliwości, w którym moc maleje proporcjonalnie do (1/f)^β, gdzie f jest częstotliwością, a β jest znana jako składowa widmowa ^47,48. Model fBm prowadzi do obrazów fraktalnych, które są podobne do losowych wzorców mętnych; Podobne wzorce są często obserwowane w mammografii. Aby zastosować ten model do obrazów, dwuwymiarowe (2D) widmo Fouriera obrazu musi zostać przekształcone w funkcję jednowymiarową (1D).

Opracowaliśmy zintegrowaną metodę do charakteryzowania rozpiętości kątowej i do oszacowania wymiaru fraktalnego (FD) obrazu, odwzorowując widmo Fouriera 2D obrazu we współrzędnych prostokątnych, oznaczonych przez S(u,v), na widmo we współrzędnych biegunowych, oznaczone przez S(ƒ,Ν). Procedura jest opisana w następujących krokach ²¹:

1. Zastosuj okno von Hanna do każdego automatycznie wyodrębnionego zwrotu z inwestycji 128 x 128 i uzupełnij wynik zerami do tablicy o rozmiarze 256 x 256 pikseli.
2. Oblicz transformację Fouriera 2D wyściełanego ROI i wielkość każdej wynikowej wartości zespolonej, aby uzyskać oszacowanie widma mocy, S(u,v), ROI.
3. Zidentyfikuj wybrane części widma o niskiej i wysokiej częstotliwości do wykluczenia w kolejnych krokach.
4. Odwzorować widmo mocy 2D S(u,v) od współrzędnych kartezjańskich (prostokątnych) (u,v) do współrzędnych biegunowych (ƒ,Ν), aby uzyskać S(ƒ,Ν), poprzez ponowne próbkowanie i obliczenie średniej ważonej czterech sąsiadów każdego punktu dla odległości radialnej f w zakresie od zera do połowy częstotliwości próbkowania oraz w zakresie kąta Ν = [0, 179°].
5. Przekształć widmo 2D S(ƒ,Ν) w funkcję 1D S(ƒ), całkując jako funkcję odległości radialnej lub częstotliwości f od punktu o zerowej częstotliwości w zakresie w Ν = [0, 179°] w kącie.
6. Zastosować regresję liniową do ograniczonego zakresu częstotliwości widma 1D S(ƒ) w skali logarytmicznej, z wyłączeniem punktów w wybranych obszarach niskiej i wysokiej częstotliwości, oraz uzyskać β nachylenia dopasowanej linii, która jest oszacowaniem składowej widmowej w modelu fBm.
7. Oblicz szacowaną wartość FD jako ^15,49,50 FD = (8 - β) / 2.
8. Przekształć widmo 2D S(ƒ,Ν) w funkcję 1D S(Ν), całkując jako funkcję kąta Ν dla zakresu [0, 179°], od punktu o zerowej częstotliwości przez odległość radialną ƒ = [1, 128] pikseli.
9. Znormalizuj S(Ν), aby uzyskać sumę jednostkową i oblicz entropię wyniku jako .

Transformacja geometryczna opisana powyżej prowadzi do lepszego przedstawienia i wizualizacji charakterystyk spektralnych okresowej lub spikulowanej tekstury ⁹. Wybrane obszary o niskiej i wysokiej częstotliwości muszą zostać wykluczone, aby usunąć wpływ składowych o niskiej częstotliwości związanych z ogólnym wyglądem obrazu i dużymi strukturami obecnymi na obrazie, a także aby zapobiec skutkom szumu o wysokiej częstotliwości. W niniejszej pracy zakresy częstotliwości, które mają być wyłączone z szacowania β i FD (tj. części nieliniowe), zostały wybrane na podstawie eksperymentów z wykorzystaniem zsyntetyzowanych obrazów ze znanymi FD, a także przy użyciu pewnej liczby ROI mammografii. Zakres ƒ używany do dopasowania modelu liniowego odpowiada [6, 96] pikselom lub [0,117, 1,875] ^mm-1, gdzie zakres [1,128] pikseli odpowiada dyskretnej reprezentacji zakresu częstotliwości [0, 2,5] ^mm-1.

Rysunki 9 i 10 ilustrują różne etapy analizy fraktalnej i szacowania kątowego rozkładu mocy w dziedzinie częstotliwości odpowiednio dla TP ROI i FP ROI. Rysunek 9D wskazuje na istnienie wielokierunkowych wzorców spikulacyjnych dla TP ROI, podczas gdy na rysunku 10D rozrzut mocy jest ograniczony do niewielkiej liczby pasm kątowych dla FP ROI.

7. Klasyfikacja i walidacja wzorca

Mamy teraz trzy miary lub cechy dla każdego automatycznie wykrytego ROI: [wartość węzła, H_F, FD], jak opisano w sekcji 6. Używamy tych cech, indywidualnie i zbiorowo, do scharakteryzowania narastających wzorców związanych ze zniekształceniami architektonicznymi oraz do odróżnienia wykrytych ROI TP od wykrytych ROI FP.

Dla zwrotu z inwestycji TP pokazanego na rysunku 9, wektor cech złożony z trzech otrzymanych miar to [0,0299, 7,2224, 2,3037]. Dla zwrotu z inwestycji w FP pokazanego na rysunku 10 odpowiadający mu wektor cech to [0,0349, 6,9444, 2,5223]. Zgodnie z oczekiwaniami wartość H_F jest wyższa dla zwrotu z inwestycji TP niż dla zwrotu z inwestycji w FP, a wartość FD jest niższa. Jednak wartość węzła jest niższa dla TP ROI niż dla FP ROI, co jest sprzeczne z oczekiwanymi różnicami, ze względu na obecność nakładających się struktur w tym ostatnim. Ogólnie rzecz biorąc, możemy oczekiwać, że niektóre funkcje będą podążać za oczekiwanymi trendami i pomogą w klasyfikacji ROI, nawet jeśli inne funkcje nie wykazują oczekiwanego zachowania.

Aby ocenić wydajność funkcji, używamy obszaru pod krzywą charakterystyki operacyjnej odbiornika (ROC) (AUC) ^51,52 i analizy ROC (FROC) ^53-55. Do analizy ROC i FROC z indywidualną cechą nie używamy wyszkolonego klasyfikatora; Zamiast tego stosujemy próg przesuwny. Trzy indywidualne cechy węzła, FD i H_F zapewniły wartości AUC odpowiednio 0,61, 0,59 i 0,64, co wskazuje na dobry potencjał, ale nieodpowiednią wydajność w klasyfikacji wzorców. Wartości p tych samych cech to 1,7638e-009, 1,8793e-004 i 2,2615e-013, co wskazuje na statystycznie bardzo istotne różnice między ich wartościami dla zestawów automatycznie wykrywanych ROI TP i FP.

Gdy do reprezentacji próbek do klasyfikacji używa się dużej liczby cech, konieczne jest wybranie optymalnego podzbioru cech, aby usunąć skorelowane cechy i zmniejszyć złożoność klasyfikatora ^20-22; kilka procedur, takich jak krokowa regresja logistyczna ⁵⁶, może być użytych do tego celu. W niniejszej pracy, ponieważ używamy tylko trzech funkcji na zwrot z inwestycji, nie dokonujemy wyboru funkcji.

Aby przeprowadzić walidację wyszkolonego klasyfikatora, stosujemy procedury klasyfikacji wzorców z podejściem "pomiń jednego pacjenta". Wykluczamy wszystkie ROI uzyskane z mammografii badanej pacjentki z procedury szkoleniowej klasyfikatora, a następnie stosujemy uzyskany w ten sposób klasyfikator do przypadku testowego. Następnie powtarzamy procedurę dla całego zbioru danych, jednego przypadku lub pacjenta na raz.

Do analizy ROC z zestawem trzech cech, używamy klasyfikatora, który przeprowadza kwadratową analizę dyskryminacyjną przy założeniu Bayesa ⁵⁷. Aby wygenerować krzywe FROC, bierzemy pod uwagę zwrot z inwestycji TP o najwyższej wartości dyskryminacyjnej w dwóch obrazach mammograficznych dostępnych dla pacjentki, z wyjątkiem sześciu przypadków, w których dostępny jest tylko jeden obraz w każdym przypadku.

Wyniki

Trzy cechy, a mianowicie wartość węzła, FD i H_F, dostarczały wartości AUC odpowiednio 0,61, 0,59 i 0,64, gdy każda funkcja była używana osobno. Połączone użycie tych trzech funkcji zapewniło lepszą wydajność przy AUC = 0,70. Krzywa FROC uzyskana w wyniku połączenia tych trzech cech pokazano na rysunku 11, który wskazuje czułość 80% przy 5,6 FPs/pacjenta i 89% przy 7,5 FPs/pacjenta. Użycie tylko wartości węzła zapewniło czułość 80% przy 8,1 FPs/pacjenta i 89% przy 13,8 FPs/pacjenta.

Dyskusja

Przedstawiliśmy szereg zaawansowanych technik cyfrowego przetwarzania obrazu i rozpoznawania wzorców, znanych również jako uczenie maszynowe i CAD, do wykrywania zniekształceń architektonicznych we wcześniejszych mammografiach przypadków raka interwałowego. Metody te opierają się na analizie zorientowanych wzorców teksturalnych obecnych w obrazach mammograficznych. Nasze metody, w tym kilka innych funkcji zaproponowanych w naszych powiązanych pracach, są w stanie wykryć wczesne objawy raka piersi na 15 miesięcy przed cza...