Recherche de régions / motifs symétriques dans l'image

14

J'ai un ensemble d'images représentant la courbure moyenne d'une surface arrière humaine.

Ce que je veux faire est de "scanner" l'image pour des points qui ont des "homologues" réfléchis similaires dans une autre partie de l'image (très probablement symétrique à la ligne médiane, mais pas nécessairement car il pourrait y avoir des déformations). Certaines techniques d'assemblage d'images l'utilisent pour «détecter automatiquement» des points similaires entre les images, mais je veux les détecter des deux côtés de la même image.

Le but ultime est de trouver une ligne longitudinale continue, très probablement incurvée, qui divise le dos de manière adaptative en «moitiés» symétriques.

Un exemple d'image est placé ci-dessous. Notez que toutes les régions ne sont pas symétriques (spécifiquement, juste au-dessus du centre de l'image, la "bande" verticale rouge dévie vers la droite). Cette région devrait recevoir un mauvais score, ou autre, mais alors la symétrie locale serait définie à partir de points symétriques placés plus loin. Dans tous les cas, je devrai adapter n'importe quel algorithme à mon domaine d'application, mais ce que je recherche c'est une stratégie de corrélation / convolution / correspondance de modèle, il doit déjà y avoir quelque chose autour, je pense.

(EDIT: il y a plus d'images ci-dessous, et quelques explications supplémentaires)

entrez la description de l'image ici

EDIT: comme demandé, je vais inclure des images plus typiques, soit bien comportées et problématiques. Mais au lieu d'images colorées, ce sont des images en niveaux de gris, de sorte que la couleur se rapporte directement à la magnitude des données, ce qui ne s'est pas produit avec l'image colorée (fournie uniquement pour la communication). Bien que les images grises semblent manquer de contraste par rapport à la couleur, les dégradés de données sont là et peuvent être affichés avec un contraste adaptatif si vous le souhaitez.


1) Image d'un sujet très symétrique:

entrez la description de l'image ici


2) Image du même sujet à un moment différent. Bien qu'il y ait plus de "fonctionnalités" (plus de dégradés), il ne se "sent" pas aussi symétrique qu'auparavant:

entrez la description de l'image ici


3) Un jeune sujet mince, avec des convexités (protubérances osseuses, désignées par des régions plus claires) à la ligne médiane au lieu de la ligne médiane concave la plus courante:

entrez la description de l'image ici


4) Un jeune avec une déviation vertébrale confirmée par radiographie (notez les asymétries):

entrez la description de l'image ici


5) Le sujet "incliné" typique (bien que généralement symétrique autour de la ligne médiane incurvée, et en tant que tel pas correctement "déformé"):

entrez la description de l'image ici


Toute aide est la bienvenue!

heltonbiker
la source
Pourquoi ne pas simplement utiliser la colonne vertébrale comme séparateur?
Jim Clay
@ JimClay: Je soupçonne que la colonne vertébrale est la partie mesurée, par rapport à l'axe de symétrie réel du reste de l'image
endolith
"Certaines techniques d'assemblage d'images l'utilisent pour" détecter automatiquement "des points similaires entre les images." Faites une copie inversée de l'image, puis utilisez-en une. :)
endolith
Ne pourriez-vous pas simplement refléter l'image le long de l'axe Y et utiliser un algorithme d'enregistrement? Parce qu'il y a déjà beaucoup de recherches sur les algorithmes d'enregistrement flexibles / non paramétriques sur lesquels vous pouvez vous appuyer.
Niki Estner
JimClay, la colonne vertébrale est ce que je veux trouver, je ne sais pas où elle est; Endolith, ma question implique que les gens me disent les noms de certains de ces algorithmes, je n'en ai pas encore trouvé. Et Nikie, c'est tout le problème, mais je ne connais aucun de ces algorithmes, c'est pourquoi je pose la question en premier lieu: o)
heltonbiker

Réponses:

9

Comme je l'ai dit dans les commentaires, l'enregistrement d'images médicales est un sujet avec beaucoup de recherches disponibles, et je ne suis pas un expert. D'après ce que j'ai lu, l'idée de base couramment utilisée est de définir un mappage entre deux images (dans votre cas, une image et son image miroir), puis de définir des termes énergétiques pour la fluidité et la similitude de l'image si le mappage est appliqué, et enfin optimiser ce mappage à l'aide de techniques d'optimisation standard (ou parfois spécifiques à l'application).

J'ai piraté un algorithme rapide dans Mathematica pour le démontrer. Ce n'est pas un algorithme que vous devez utiliser dans une application médicale, seulement une démonstration des idées de base.

Tout d'abord, je charge votre image, la reflète et divise ces images en petits blocs:

src = ColorConvert[Import["http://i.stack.imgur.com/jf709.jpg"], 
   "Grayscale"];
mirror = ImageReflect[src, Left -> Right];
blockSize = 30;
partsS = ImagePartition[src, {blockSize, blockSize}];
partsM = ImagePartition[mirror, {blockSize, blockSize}];
GraphicsGrid[partsS]

Graphiques Mathematica

Normalement, nous ferions un enregistrement rigide approximatif (en utilisant par exemple des points clés ou des moments d'image), mais votre image est presque centrée, donc je vais sauter ceci.

Si nous regardons un bloc et son homologue d'image miroir:

{partsS[[6, 10]], partsM[[6, 10]]}

Graphiques Mathematica

Nous pouvons voir qu'ils sont similaires, mais décalés. La quantité et la direction du changement sont ce que nous essayons de découvrir.

Pour quantifier la similitude de correspondance, je peux utiliser la distance euclidienne au carré:

ListPlot3D[
  ImageData[
   ImageCorrelate[partsM[[6, 10]], partsS[[6, 10]], 
    SquaredEuclideanDistance]]]

Graphiques Mathematica

malheureusement, l'utilisation de ces données est l'optimisation directement plus difficile que je ne le pensais, j'ai donc utilisé une approximation de second ordre à la place:

fitTerms = {1, x, x^2, y, y^2, x*y};

fit = Fit[
   Flatten[MapIndexed[{#2[[1]] - blockSize/2, #2[[2]] - 
        blockSize/2, #1} &, 
     ImageData[
      ImageCorrelate[partsM[[6, 10]], partsS[[6, 10]], 
       SquaredEuclideanDistance]], {2}], 1], fitTerms, {x, y}];

Plot3D[fit, {x, -25, 25}, {y, -25, 25}]

Graphiques Mathematica

La fonction n'est pas la même que la fonction de corrélation réelle, mais elle est suffisamment proche pour une première étape. Calculons ceci pour chaque paire de blocs:

distancesFit = MapThread[
   Function[{part, template},
    Fit[Flatten[
      MapIndexed[{#2[[2]] - blockSize/2, #2[[1]] - blockSize/2, #1} &,
        ImageData[
        ImageCorrelate[part, template, 
         SquaredEuclideanDistance]], {2}], 1], 
     fitTerms, {x, y}]], {partsM, partsS}, 2];

Cela nous donne notre premier terme énergétique pour l'optimisation:

variablesX = Array[dx, Dimensions[partsS]];
variablesY = Array[dy, Dimensions[partsS]];

matchEnergyFit = 
  Total[MapThread[#1 /. {x -> #2, y -> #3} &, {distancesFit, 
     variablesX, variablesY}, 2], 3];

variablesX/Ycontient les décalages pour chaque bloc et donne une matchEnergyFitapproximation de la différence euclidienne au carré entre l'image d'origine et l'image miroir avec les décalages appliqués.

L'optimisation de cette énergie seule donnerait de mauvais résultats (si elle convergeait du tout). Nous voulons également que les décalages soient lisses, où la similitude des blocs ne dit rien sur le décalage (par exemple le long d'une ligne droite ou sur le fond blanc).

Nous avons donc mis en place un deuxième terme énergétique pour la douceur:

smoothnessEnergy = Total[Flatten[
    {
     Table[
      variablesX[[i, j - 1]] - 2 variablesX[[i, j]] + 
       variablesX[[i, j + 1]], {i, 1, Length[partsS]}, {j, 2, 
       Length[partsS[[1]]] - 1}],
     Table[
      variablesX[[i - 1, j]] - 2 variablesX[[i, j]] + 
       variablesX[[i + 1, j]], {i, 2, Length[partsS] - 1}, {j, 1, 
       Length[partsS[[1]]]}],
     Table[
      variablesY[[i, j - 1]] - 2 variablesY[[i, j]] + 
       variablesY[[i, j + 1]], {i, 1, Length[partsS]}, {j, 2, 
       Length[partsS[[1]]] - 1}],
     Table[
      variablesY[[i - 1, j]] - 2 variablesY[[i, j]] + 
       variablesY[[i + 1, j]], {i, 2, Length[partsS] - 1}, {j, 1, 
       Length[partsS[[1]]]}]
     }^2]];

Heureusement, l'optimisation contrainte est intégrée dans Mathematica:

allVariables = Flatten[{variablesX, variablesY}];
constraints = -blockSize/3. < # < blockSize/3. & /@ allVariables;
initialValues = {#, 0} & /@ allVariables;
solution = 
  FindMinimum[{matchEnergyFit + 0.1 smoothnessEnergy, constraints}, 
   initialValues];

Regardons le résultat:

grid = Table[{(j - 0.5)*blockSize - dx[i, j], (i - 0.5)*blockSize - 
      dy[i, j]}, {i, Length[partsS]}, {j, Length[partsS[[1]]]}] /. 
   solution[[2]];
Show[src, Graphics[
  {Red,
   Line /@ grid,
   Line /@ Transpose[grid]
   }]]

Graphiques Mathematica

Le 0.1facteur avant smoothnessEnergyest le poids relatif que l'énergie de lissage obtient par rapport au terme d'énergie de correspondance d'image. Ce sont des résultats pour différents poids:

Graphiques Mathematica

Améliorations possibles:

  • Comme je l'ai dit, effectuez d'abord un enregistrement rigide. Avec un fond blanc, un enregistrement basé sur des moments d'image simples devrait fonctionner correctement.
  • Ce n'est qu'une étape. Vous pouvez utiliser les décalages que vous avez trouvés dans une étape et les améliorer dans une deuxième étape, peut-être avec une fenêtre de recherche plus petite ou des tailles de bloc plus petites
  • J'ai lu des articles où ils le font sans blocs du tout, mais optimisent un décalage par pixel.
  • Essayez différentes fonctions de douceur
Niki Estner
la source
répondre trop longtemps pour lire juste pour le plaisir, mais l'image finale est assez indicative: elle a l'air incroyable: D
penelope
Cette réponse a été très éclairante. J'aurai besoin de temps pour l'avaler, mais très probablement la technique d'enregistrement non rigide est ce que je devrai utiliser. Heureusement, vous avez fourni quelques détails conceptuels, donc dans le pire des cas, je peux trouver une approche similaire. En attendant, je mettrai à jour la question avec plus d'images. Merci pour l'instant!
heltonbiker
4

Question interessante. Tout d'abord, vous recherchez peut-être des approches basées sur le détecteur de points d'intérêt et la correspondance. Cela comprendrait SIFT (Scale-Invariant Feature Transform), SURF, ORB, etc ... ou même une approche plus simple basée uniquement sur l'opérateur Harris (csce.uark.edu/~jgauch/library/Features/Harris.1988.pdf ). Il n'est pas clair de votre message ce que vous avez essayé, donc je suis désolé si je suis naïf ici.

Dit cela, permettez-moi d'adopter une approche plus simple avec la morphologie mathématique (MM) juste pour le plaisir :) Les images pour la visualisation de toutes les étapes sont à la fin.

J'ai pris votre exemple d'image et l'ai converti en espace colorimétrique L a b * à l'aide d'ImageMagick et j'ai utilisé uniquement la bande L *:

convert x.jpg -colorspace Lab -separate %d.png

0.png correspond à la bande L *. Maintenant, je suis sûr que vous avez les données d'image réelles, mais je traite des artefacts de compression jpg et ce qui ne l'est pas. Pour gérer partiellement ce problème, j'ai effectué une ouverture morphologique suivie d'une fermeture morphologique avec un disque plat de rayon 5. C'est un moyen de base pour réduire le bruit avec MM, et étant donné le rayon du disque, peu d'image est modifiée. Ensuite, mon idée était basée sur cette image unique, qui a de grandes chances d'échouer pour d'autres cas. Votre région d'intérêt se distingue visuellement en étant plus sombre ("plus chaude" dans votre image couleur), j'ai donc supposé qu'un binariseur basé sur des statistiques pourrait bien fonctionner. J'ai utilisé l'approche d'Otsu, qui est automatique.

À ce stade, il est possible de visualiser clairement la région centrale d'intérêt. Le problème est que, dans mon approche, je voulais que ce soit un composant fermé, mais ce n'est pas le cas. Je commence par éliminer tous les composants connectés qui sont plus petits que le plus grand (sans compter l'arrière-plan comme l'un d'eux). Cela a plus de chance de fonctionner dans d'autres cas si le résultat de la binarisation était bon. Dans votre exemple d'image, un composant est connecté à l'arrière-plan, il n'est donc pas ignoré, mais cela ne pose pas de problème.

Si vous me suivez toujours, nous n'avons pas encore trouvé la véritable région centrale supposée d'intérêt. Voici mon point de vue à ce sujet. Peu importe la courbe de la personne (en fait, je peux voir certains cas problématiques), la région ressemble à une ligne verticale. À cette fin, je simplifie l'image actuelle en effectuant une ouverture morphologique avec une ligne verticale de longueur 100. Cette longueur est purement arbitraire, si vous n'avez pas de problèmes de mise à l'échelle, ce n'est pas une valeur difficile à déterminer. Maintenant, nous rejetons à nouveau les composants, mais j'ai été un peu plus prudent à cette étape. J'ai utilisé l'ouverture par zone avec le complément de l'image pour éliminer ce que je considérais comme de petites régions, cela pourrait être fait de manière plus contrôlée en effectuant quelque chose sous la forme d'une analyse granulométrique (également de MM).

Nous avons à peu près trois morceaux maintenant: la partie gauche de l'image, la partie centrale et la partie droite de l'image. On s'attend à ce que la partie centrale soit la plus petite composante des trois, de sorte qu'elle est obtenue de manière triviale.

Voici le résultat final, l'image en bas à droite est juste l'image superposée à sa gauche avec l'original. Les chiffres individuels ne sont pas tous alignés, désolé pour la hâte.

http://i.imgur.com/XRhYv.png

mmgp
la source
Merci beaucoup pour votre aimable intérêt, mais votre approche devrait prendre en compte certaines propriétés de mes données (pas une plainte, juste un détail): 1) Les données réelles sont un tableau 2D de flotteurs, assortis d'une couleur rouge-jaune divergente- palette de couleurs verte dans le matplotlib de Python. Je ne pense pas que travailler avec les données de couleur serait conceptuellement correct, les images sont montrées uniquement à des fins de communication; 2) Les données réelles concernent la courbure de surface (convexe contre concave), les parties rouges étant concaves et les vertes étant convexes. L'axe symétrique ne tombe pas nécessairement dans une région concave.
heltonbiker
Je vais très prochainement ajouter quelques images (et remplacer celle-ci) avec des images en niveaux de gris, afin que les images elles-mêmes puissent être utilisées pour les tests, éliminant ainsi le risque de distorsion de la plage dynamique due à la couleur.
heltonbiker
Malheureusement, les données ne sont pas encore disponibles. Les images en niveaux de gris en sont au mieux une approximation.
mmgp
Je pense que l'approximation est probablement suffisante, mais cela ne me dérange pas de fournir les données réelles. Je peux publier des liens de téléchargement publics DropBox, mais je ne sais pas dans quel format de fichier.
heltonbiker