1

Tengo un conjunto de datos que he extraído de un DataFrame hecho con pandas, formado por 18383 registros (x,y,z0,z,arbol). Dichos datos se corresponden con la sección de un árbol y por eso se distribuyen en el espacio en forma de círculo. Me gustaría hacer un ajuste por mínimos cuadrados para obtener el diámetro del círculo que describen. Soy muy novata en Python y encontrado este fragmento de código que permite hacer el ajuste, sin embargo no sé cómo introducir mi DataFrame en él para hacer el ajuste, ¿debería convertirlo en una matriz?. El código es:

method_2 = "leastsq"
from numpy import *
from scipy import optimize

def calc_R(xc, yc):
    """ calculate the distance of each 2D points from the center (xc, yc) """
    return sqrt((x-xc)**2 + (y-yc)**2)

def f_2(c):
    """ calculate the algebraic distance between the data points and the mean circle centered at c=(xc, yc) """
    Ri = calc_R(*c)
    return Ri - Ri.mean()

center_estimate = x_m, y_m
center_2, ier = optimize.leastsq(f_2, center_estimate)

xc_2, yc_2 = center_2
Ri_2       = calc_R(*center_2)
R_2        = Ri_2.mean()
residu_2   = sum((Ri_2 - R_2)**2)

Adjunto una imagen del gráfico que he hecho con mis datos así como de los datos que estoy utilizando

8
  • podrias compartir el dataframe como texto, no pensaras que vamos escribir esos numeros
    – user22539
    Commented el 17 jul. 2018 a las 8:27
  • Perdona, todavía no manejo muy bien la página...Ahí va el enlace de los datos en formato .txt dropbox.com/s/ut3mi8etuqmo7wh/0arbol1.7z?dl=0
    – Deva
    Commented el 17 jul. 2018 a las 8:53
  • entonces lee Cómo preguntar por favor y pasa el recorrido de bienvenida
    – user22539
    Commented el 17 jul. 2018 a las 8:54
  • @Covadonga El enlace que has puesto lleva a un archivo .7z que está vacío.
    – abulafia
    Commented el 17 jul. 2018 a las 9:46
  • @Covadonga Por otro lado, tu pregunta se parece mucho a esta otra que respondí hace un tiempo. Mira si te puede ser de utilidad mi respuesta.
    – abulafia
    Commented el 17 jul. 2018 a las 9:48

1 respuesta 1

0

Basándome en mi respuesta a otra pregunta sobre ajuste de círculos, la adaptación es muy sencilla.

En lugar las listas x e y que usaba el problema original para contener las coordenadas de los puntos, en esta ocasión usaremos arbol.x y arbol.y, siendo arbol un dataframe pandas que contiene los datos de un árbol (y que tiene al menos dos columnas llamadas x e y).

El siguiente código muestra cómo operar (he ejecutado todos los comandos mostrados en celdas de un Jupyter Notebook):

  1. Bajarse los datos proporcionados por la usuaria

    !wget "https://www.dropbox.com/s/mc7njkif3ql78ts/arbol1.txt?dl=1" -o /dev/null -O arbol1.txt
    
  2. Leerlos en un dataframe pandas y ver que tienen buen aspecto

    # Leer los datos con pandas
    import pandas as pd
    df = pd.read_table("arbol1.txt", header=0)
    df.head()
    

            x       y     z0       z  arbol
    0  23.215  73.319  1.210  1068.3      1
    1  23.209  73.332  1.210  1068.3      1
    2  23.203  73.344  1.209  1068.3      1
    3  23.208  73.321  1.397  1068.5      1
    4  23.207  73.321  1.369  1068.4      1
    
  3. Quedarse con los datos del "arbol 1" y pintarlos para ver que siguen una forma más o menos circular:

    # Quedarse con un arbol, por ejemplo el 3
    arbol = df[df.arbol == 1]
    ax = arbol.plot(x="x", y="y", kind="scatter")
    ax.set_aspect('equal')
    

    introducir la descripción de la imagen aquí

  4. Utilizar el código de la otra pregunta, adaptado a este caso, para calcular el centro y radio del círculo:

    # Hacer el ajuste (basado en otra respuesta mía)
    # https://es.stackoverflow.com/a/166354/7123
    
    from scipy import optimize as opt
    import numpy as np
    
    def calc_R(x, y, xc, yc):
        """ Calculate the distance of each 2D points from the center (xc, yc) """
        return np.sqrt((x-xc)**2 + (y-yc)**2)
    
    def f_2(c, x, y):
        """ Calculate the algebraic distance between the 2D points and the mean circle centered at c=(xc, yc) """
        Ri = calc_R(x, y, *c)
        return Ri - Ri.mean()
    
    # Coordenadas del baricentro (estimación inicial del centro del círculo)
    x_m = arbol.x.mean()
    y_m = arbol.y.mean()
    barycenter = x_m, y_m
    # Aplicamos un algoritmo para aproximar el centro de la circunferencia.
    circle_center, ier = opt.leastsq(f_2, barycenter, args=(arbol.x, arbol.y))
    xc, yc = circle_center
    print("Centro del círculo:", xc, yc)
    
    # El radio encontrado
    r = calc_R(arbol.x, arbol.y, *circle_center).mean()
    print("Radio del árbol:", r)
    

    Centro del círculo: 22.98550108357253 73.229096837513
    Radio del árbol: 0.2332831108455863
    

Finalmente podemos pintar encima de la gráfica anterior el círculo que hemos encontrado. Así:

from matplotlib.patches import Circle
from matplotlib.collections import PatchCollection

ax = arbol.plot(x="x", y="y", kind="scatter")
ax.plot(x=xc, y=yc, mark="*", color="red")
ax.add_collection(PatchCollection([Circle((xc, yc), r)], color="orange", alpha=0.2))
ax.set_aspect('equal')
ax.set_xlim((xc-r*1.1, xc+r*1.1))
ax.set_ylim((yc-r*1.1, yc+r*1.1))

Resultado

Tu Respuesta

By clicking “Publica tu respuesta”, you agree to our terms of service and acknowledge you have read our privacy policy.

¿No es la respuesta que buscas? Examina otras preguntas con la etiqueta o formula tu propia pregunta.