Détection des pics dans un réseau 2D

Question

J'aide une clinique vétérinaire à mesurer la pression sous la patte d'un chien. J'utilise Python pour l'analyse de mes données et je suis maintenant coincé en essayant de diviser les pattes en sous-régions (anatomiques).

J'ai fait un tableau 2D de chaque patte, qui consiste en les valeurs maximales pour chaque capteur qui a été chargé par la patte au fil du temps. Voici un exemple d'une patte, où j'ai utilisé Excel pour dessiner les zones que je veux "détecter". Il s'agit de cases 2 par 2 autour du capteur avec des maxima locaux, qui ensemble ont la plus grande somme.

J'ai donc fait quelques expériences et j'ai décidé de chercher simplement les maximums de chaque colonne et ligne (je ne peux pas regarder dans une seule direction à cause de la forme de la patte). Cela semble " détecter " assez bien l'emplacement des orteils séparés, mais cela marque aussi les capteurs voisins.

Quelle serait donc la meilleure façon d'indiquer à Python quels sont les maximums que je souhaite obtenir ?

Remarque : les carrés 2x2 ne peuvent pas se chevaucher, car ils doivent être des orteils séparés !

J'ai également pris 2x2 par commodité, toute solution plus avancée est la bienvenue, mais je suis simplement un scientifique du mouvement humain, donc je ne suis ni un vrai programmeur ni un mathématicien, alors s'il vous plaît, restez "simple".

Voici un qui peut être chargé avec np.loadtxt

---

Résultats

J'ai donc essayé la solution de @jextee (voir les résultats ci-dessous). Comme vous pouvez le voir, ça marche très bien sur les pattes avant, mais ça marche moins bien pour les pattes arrière.

Plus précisément, il ne peut pas reconnaître le petit pic qui constitue le quatrième orteil. Ceci est évidemment inhérent au fait que la boucle regarde de haut en bas vers la valeur la plus basse, sans tenir compte de l'endroit où elle se trouve.

Quelqu'un saurait-il comment modifier l'algorithme de @jextee pour qu'il soit capable de trouver le quatrième orteil ?

Comme je n'ai pas encore traité d'autres essais, je ne peux pas fournir d'autres échantillons. Mais les données que j'ai données avant étaient les moyennes de chaque patte. Ce fichier est un tableau avec les données maximales de 9 pattes dans l'ordre où elles sont entrées en contact avec la plaque.

Cette image montre comment ils ont été spatialement répartis sur la plaque.

Mise à jour :

J'ai créé un blog pour toute personne intéressée y J'ai créé un OneDrive avec toutes les mesures brutes. Donc, à tous ceux qui demandent plus de données : plus de pouvoir pour vous !

---

Nouvelle mise à jour :

Donc, après l'aide que j'ai reçue pour mes questions concernant détection des pattes y triage des pattes J'ai enfin pu vérifier la détection des orteils pour chaque patte ! Il s'avère que ça ne fonctionne pas si bien que ça, sauf pour les pattes de la taille de celle de mon exemple. Bien sûr, avec le recul, c'est ma propre faute pour avoir choisi le 2x2 de façon si arbitraire.

Voici un bel exemple de ce qui ne va pas : un ongle est reconnu comme un orteil et le " talon " est si large qu'il est reconnu deux fois !

La patte est trop grande, donc en prenant une taille 2x2 sans chevauchement, certains orteils sont détectés deux fois. Dans l'autre sens, chez les petits chiens, il arrive souvent que le cinquième orteil ne soit pas détecté, ce que je soupçonne d'être dû au fait que la zone 2x2 est trop grande.

Après j'essaie la solution actuelle sur toutes mes mesures Je suis arrivé à la conclusion stupéfiante que pour presque tous mes petits chiens, il n'a pas trouvé de 5e orteil et que dans plus de 50 % des impacts pour les grands chiens, il en trouverait davantage !

Il est donc clair que je dois le changer. J'ai pensé qu'il fallait changer la taille de la neighborhood à quelque chose de plus petit pour les petits chiens et de plus grand pour les grands chiens. Mais generate_binary_structure ne me laissait pas changer la taille du tableau.

Par conséquent, j'espère que quelqu'un d'autre a une meilleure suggestion pour situer les orteils, peut-être en faisant en sorte que la surface des orteils soit proportionnelle à la taille de la patte ?

Answer 1

J'ai détecté les pics en utilisant un filtre local maximal . Voici le résultat sur votre premier jeu de données de 4 pattes :

Je l'ai également exécuté sur le deuxième ensemble de données de 9 pattes et ça a marché aussi bien .

Voici comment procéder :

import numpy as np
from scipy.ndimage.filters import maximum_filter
from scipy.ndimage.morphology import generate_binary_structure, binary_erosion
import matplotlib.pyplot as pp

#for some reason I had to reshape. Numpy ignored the shape header.
paws_data = np.loadtxt("paws.txt").reshape(4,11,14)

#getting a list of images
paws = [p.squeeze() for p in np.vsplit(paws_data,4)]

def detect_peaks(image):
    """
    Takes an image and detect the peaks usingthe local maximum filter.
    Returns a boolean mask of the peaks (i.e. 1 when
    the pixel's value is the neighborhood maximum, 0 otherwise)
    """

    # define an 8-connected neighborhood
    neighborhood = generate_binary_structure(2,2)

    #apply the local maximum filter; all pixel of maximal value 
    #in their neighborhood are set to 1
    local_max = maximum_filter(image, footprint=neighborhood)==image
    #local_max is a mask that contains the peaks we are 
    #looking for, but also the background.
    #In order to isolate the peaks we must remove the background from the mask.

    #we create the mask of the background
    background = (image==0)

    #a little technicality: we must erode the background in order to 
    #successfully subtract it form local_max, otherwise a line will 
    #appear along the background border (artifact of the local maximum filter)
    eroded_background = binary_erosion(background, structure=neighborhood, border_value=1)

    #we obtain the final mask, containing only peaks, 
    #by removing the background from the local_max mask (xor operation)
    detected_peaks = local_max ^ eroded_background

    return detected_peaks

#applying the detection and plotting results
for i, paw in enumerate(paws):
    detected_peaks = detect_peaks(paw)
    pp.subplot(4,2,(2*i+1))
    pp.imshow(paw)
    pp.subplot(4,2,(2*i+2) )
    pp.imshow(detected_peaks)

pp.show()

Tout ce que vous devez faire après est d'utiliser scipy.ndimage.measurements.label sur le masque pour étiqueter tous les objets distincts. Vous pourrez ensuite jouer avec eux individuellement.

Note que la méthode fonctionne bien parce que le fond n'est pas bruyant. S'il l'était, vous détecteriez un tas d'autres pics indésirables dans l'arrière-plan. Un autre facteur important est la taille de la quartier . Vous devrez l'ajuster si la taille du pic change (elle doit rester à peu près proportionnelle).

Answer 2

Solution

Fichier de données : paw.txt . Code source :

from scipy import *
from operator import itemgetter

n = 5  # how many fingers are we looking for

d = loadtxt("paw.txt")
width, height = d.shape

# Create an array where every element is a sum of 2x2 squares.

fourSums = d[:-1,:-1] + d[1:,:-1] + d[1:,1:] + d[:-1,1:]

# Find positions of the fingers.

# Pair each sum with its position number (from 0 to width*height-1),

pairs = zip(arange(width*height), fourSums.flatten())

# Sort by descending sum value, filter overlapping squares

def drop_overlapping(pairs):
    no_overlaps = []
    def does_not_overlap(p1, p2):
        i1, i2 = p1[0], p2[0]
        r1, col1 = i1 / (width-1), i1 % (width-1)
        r2, col2 = i2 / (width-1), i2 % (width-1)
        return (max(abs(r1-r2),abs(col1-col2)) >= 2)
    for p in pairs:
        if all(map(lambda prev: does_not_overlap(p,prev), no_overlaps)):
            no_overlaps.append(p)
    return no_overlaps

pairs2 = drop_overlapping(sorted(pairs, key=itemgetter(1), reverse=True))

# Take the first n with the heighest values

positions = pairs2[:n]

# Print results

print d, "\n"

for i, val in positions:
    row = i / (width-1)
    column = i % (width-1)
    print "sum = %f @ %d,%d (%d)" % (val, row, column, i)
    print d[row:row+2,column:column+2], "\n"

Sortie sans carrés superposés. Il semble que les mêmes zones soient sélectionnées que dans votre exemple.

Quelques commentaires

La partie la plus délicate consiste à calculer les sommes de tous les carrés 2x2. Je suppose que vous avez besoin de tous les carrés, donc il peut y avoir des chevauchements. J'ai utilisé des tranches pour couper les premières/dernières colonnes et lignes du tableau 2D original, puis je les ai superposées toutes ensemble et j'ai calculé les sommes.

Pour mieux comprendre, imaginez un tableau 3x3 :

>>> a = arange(9).reshape(3,3) ; a
array([[0, 1, 2],
       [3, 4, 5],
       [6, 7, 8]])

Ensuite, vous pourrez prendre ses tranches :

>>> a[:-1,:-1]
array([[0, 1],
       [3, 4]])
>>> a[1:,:-1]
array([[3, 4],
       [6, 7]])
>>> a[:-1,1:]
array([[1, 2],
       [4, 5]])
>>> a[1:,1:]
array([[4, 5],
       [7, 8]])

Imaginez maintenant que vous les empilez les uns sur les autres et que vous additionnez les éléments aux mêmes positions. Ces sommes seront exactement les mêmes sommes sur les carrés 2x2 avec le coin supérieur gauche dans la même position :

>>> sums = a[:-1,:-1] + a[1:,:-1] + a[:-1,1:] + a[1:,1:]; sums
array([[ 8, 12],
       [20, 24]])

Lorsque vous avez les sommes sur des carrés de 2x2, vous pouvez utiliser max pour trouver le maximum, ou sort ou sorted pour trouver les pics.

Pour se souvenir des positions des pics, je couple chaque valeur (la somme) avec sa position ordinale dans un tableau aplati (voir zip ). Ensuite, je calcule à nouveau la position des lignes/colonnes lorsque j'imprime les résultats.

Notes

J'ai permis aux carrés de 2x2 de se chevaucher. La version éditée filtre certains d'entre eux de sorte que seuls les carrés qui ne se chevauchent pas apparaissent dans les résultats.

Choisir ses doigts (une idée)

Un autre problème est de savoir comment choisir ce qui est susceptible d'être des doigts parmi tous les pics. J'ai une idée qui peut fonctionner ou non. Je n'ai pas le temps de la mettre en œuvre pour l'instant, donc juste un pseudo-code.

J'ai remarqué que si les doigts avant restent sur un cercle presque parfait, le doigt arrière devrait être à l'intérieur de ce cercle. En outre, les doigts avant sont plus ou moins également espacés. Nous pouvons essayer d'utiliser ces propriétés heuristiques pour détecter les doigts.

Pseudo-code :

select the top N finger candidates (not too many, 10 or 12)
consider all possible combinations of 5 out of N (use itertools.combinations)
for each combination of 5 fingers:
    for each finger out of 5:
        fit the best circle to the remaining 4
        => position of the center, radius
        check if the selected finger is inside of the circle
        check if the remaining four are evenly spread
        (for example, consider angles from the center of the circle)
        assign some cost (penalty) to this selection of 4 peaks + a rear finger
        (consider, probably weighted:
             circle fitting error,
             if the rear finger is inside,
             variance in the spreading of the front fingers,
             total intensity of 5 peaks)
choose a combination of 4 peaks + a rear peak with the lowest penalty

Il s'agit d'une approche par la force brute. Si N est relativement petit, alors je pense que c'est faisable. Pour N=12, il y a C_12^5 = 792 combinaisons, fois 5 façons de sélectionner un doigt arrière, donc 3960 cas à évaluer pour chaque patte.

Answer 3

Il s'agit d'un problème d'enregistrement des images . La stratégie générale est la suivante :

• Avoir un exemple connu, ou une sorte de avant sur les données.
• Adaptez vos données à l'exemple, ou adaptez l'exemple à vos données.
• Il est utile que vos données soient à peu près alignés en premier lieu.

Voici une approche approximative "la chose la plus stupide qui pourrait fonctionner" :

• Commencez par cinq coordonnées d'orteil à peu près à l'endroit prévu.
• Avec chacun d'entre eux, grimpez itérativement vers le sommet de la colline, c'est-à-dire, étant donné la position actuelle, déplacez-vous vers le pixel voisin maximum, si sa valeur est supérieure au pixel actuel. Arrêtez-vous lorsque les coordonnées de vos orteils ont cessé de bouger.

Pour contrer le problème d'orientation, vous pourriez disposer d'environ 8 paramètres initiaux pour les directions de base (nord, nord-est, etc.). Exécutez chacun d'eux individuellement et jetez tous les résultats où deux orteils ou plus se retrouvent au même pixel. Je vais encore y réfléchir, mais ce genre de chose fait encore l'objet de recherches dans le domaine du traitement de l'image - il n'y a pas de bonnes réponses !

Une idée un peu plus complexe : Le regroupement K-means (pondéré). Ce n'est pas si grave.

• Commencez par cinq coordonnées d'orteils, mais il s'agit maintenant de "centres de regroupement".

Puis itérer jusqu'à convergence :

• Assignez chaque pixel au cluster le plus proche (faites simplement une liste pour chaque cluster).
• Calculez le centre de masse de chaque groupe. Pour chaque cluster, c'est : Somme(coordonnée * valeur d'intensité)/Somme(coordonnée)
• Déplacez chaque groupe vers le nouveau centre de masse.

Cette méthode donnera presque certainement de bien meilleurs résultats, et vous obtiendrez la masse de chaque grappe, ce qui peut aider à identifier les orteils.

(Encore une fois, vous avez spécifié le nombre de clusters au départ. Avec le clustering, vous devez spécifier la densité d'une manière ou d'une autre : Soit vous choisissez le nombre de clusters, approprié dans ce cas, soit vous choisissez un rayon de cluster et vous voyez combien vous obtenez. Voici un exemple de cette dernière solution décalage moyen .)

Désolé pour le manque de détails de mise en œuvre ou d'autres spécificités. J'aurais bien codé tout ça, mais j'ai un délai à respecter. Si rien d'autre n'a fonctionné d'ici la semaine prochaine, faites-le moi savoir et je tenterai le coup.

Answer 4

Ce problème a été étudié de manière assez approfondie par les physiciens. Il existe une bonne implémentation dans Racine . Regardez le TSpectrum (en particulier TSpectrum2 pour votre cas) et la documentation les concernant.

Références :

1. M.Morhac et al : Méthodes d'élimination du bruit de fond pour les spectres gamma de coïncidence multidimensionnels. Nuclear Instruments and Methods in Physics Research A 401 (1997) 113-132.
2. M.Morhac et al : Efficient one- and two-dimensional Gold deconvolution and its application to gamma-ray spectra decomposition. Nuclear Instruments and Methods in Physics Research A 401 (1997) 385-408.
3. M.Morhac et al : Identification des pics dans les spectres gamma de coïncidence multidimensionnelle. Nuclear Instruments and Methods in Research Physics A 443(2000), 108-125.

...et pour ceux qui n'ont pas accès à un abonnement au NIM :

• Spectre.doc
• SpectrumDec.ps.gz
• SpectrumSrc.ps.gz
• SpectreBck.ps.gz

Answer 5

Ce ne sont que quelques idées qui me viennent à l'esprit :

• prenez le gradient (dérivée) du scan, voyez si cela élimine les faux appels.
• prendre le maximum des maxima locaux

Vous pouvez également jeter un coup d'œil à OpenCV Il est doté d'une API Python assez décente et peut proposer des fonctions qui vous seront utiles.