0

Llevo varios días peleándome con Python intentando dar con una solución para mi problema, pero me veo forzado a recurrir a la comunidad. Antes que nada decir que soy novato en Python. Os presento mi base de datos, la cual leo con pandas para obtener un dataframe. Esta dataframe se ha ordenado en función de la columna "mmsi" y la columna "time"

introducir la descripción de la imagen aquí

Tras una serie de operaciones, mi base de datos queda tal que así:

introducir la descripción de la imagen aquí

Básicamente, estas operaciones han servido para eliminar aquellos barcos (pues los datos son posiciones AIS) cuya velocidad en nudos era mayor que un nudo.

Mi propósito es utilizar estos datos para clustering, principalmente con DBSCAN.

Lo que busco es evitar que un barco solamente, que haya transmitido posiciones muy cercanas entre si durante un largo periodo de tiempo, me genere un cluster por sí mismo.

Para ello, se me había ocurrido comparar dos posiciones separadas en el tiempo (Ej: Aprox dos horas = unos 60 índices de diferencia), y ver a qué distancia se encuentran para ver si supera un umbral.

Ej: si la diferencia es mayor de 1 (millas náuticas) nos interesa esa posición, e ignoraremos las del intervalo 2-59. Es decir, nos quedaríamos con la primera, y la sesenta. El problema es que no se hacer que itere todo el dataframe de esta manera.

Acepto propuestas alternativas para evitar que un solo barco genere un cluster.

Mi intento quedó en lo siguiente:

for i in range(60,len(X)):
if X.loc[i,"mmsi"] == X.loc[i-60,"mmsi"]:
    distancia = distancia_puntos(i,i-60)
    
    if distancia < 1:
        X = X.drop(X.index[i-1:i-59])
        X = X.reset_index(drop=True)     
8
  • 1
    Un barco puede ir y volver en 2 horas. Sacaría el promedio de todas las lecturas, y si cada lectura individual está a menos de 1 milla del promedio, consideraría el barco inmóvil.
    – Candid Moe
    Commented el 3 mar. 2021 a las 17:49
  • ¿Me podrías guiar con eso? :) Commented el 3 mar. 2021 a las 18:47
  • Conozco muy poco de pandas, pero si publicas un archivo de datos manejable y representativo (con barcos que se mueven y anclados) podría intentarlo. Tendrías que ponerlo en algún sitio público y publicar el link.
    – Candid Moe
    Commented el 3 mar. 2021 a las 18:53
  • mediafire.com/file/2vqwuenw10y870g/trozo2.csv/file Aquí dejo el archivo que me pides. También cabe mencionar que estos datos han sido volcados en un archivo .csv desde una base de datos mayor. A priori, todos los barcos de este archivo están parados (el filtro fue otra columna llamada "SOG", de la cual se extrajeron los que eran igual a 0) Commented el 3 mar. 2021 a las 19:00
  • ¿mmsi es la identificación del barco?
    – Candid Moe
    Commented el 3 mar. 2021 a las 19:03

2 respuestas 2

2

Esta respuesta calcula la distancia que se mueve cada barco, sin considerar el tiempo transcurrido.

Para cada barco calcula:

  1. Posición promedio en base a su lista de coordenadas.
  2. Posición extrema, la coordenada del punto más alejado de la posición promedio.
  3. La distancia entre ambos.

Dados estos tres datos, es posible reducir su lista de coordenadas a sólo dos puntos: la posición promedio y la posición extrema.

Para calcular la distancia se ocupa la formula euclidiana ajustada según la latitud del barco.

Los resultados se pueden verificar usando calculadoras en línea como Latitude/Longitude Distance Calculator. Los resultados pueden variar, dado que el sitio ocupa la formula más precisa y general (Haversine).

No me queda claro como quieres la salida, por lo que no produce ningún archivo. Se puede fácilmente producir uno con los barcos inmóviles y sus dos coordenadas, pero no tengo claro que pasa con los que si se mueven. Espero tus comentarios.

Código

import csv
from math import sqrt, cos

miles_per_degree = 59.71  # En el ecuador

def small_distance(lat, lng, lat0, lng0):
    """
    Calcula distancia entre dos coordenadas para pequeñas distancias.
    Se aplica la formula euclidiana para distancie entre dos puntos,
    corregidas por la latitud.

    :param lat: Latitud punto partida
    :param lng: Longitud punto de partida.
    :param lat0: Latitud punto llegada
    :param lng0: Longitud punto llegada.
    :return: Distancia en millas nauticas.
    """
    x = lat - lat0
    y = (lng - lng0) * cos(lat0)
    return miles_per_degree * sqrt(x * x + y * y)


def max_distance(list_coord):
    """
    Calcula la dispersion de un conjunto de coordenadas.
    Calcula el promedio de una nube de coordenadas y luego
    la distancia de cada punto respecto al promedio.
    :param list_coord: Lista de coordenadas (lon, lat)
    :return: La maxima distancia en millas naúticas, la
             coord promedio y la coord más lejana
    """
    lon_average = sum([x[0] for x in list_coord]) / len(list_coord)
    lat_average = sum([x[1] for x in list_coord]) / len(list_coord)
    print("    coord average =", (lon_average, lat_average))
    max_dist = 0
    max_coord = None
    for (lon, lat) in list_coord:
        distance = small_distance(lat, lon, lat_average, lon_average)
        if distance >= max_dist:
            max_dist = distance
            max_coord = (lon, lat)

    print("    coord fartest =", max_coord)
    return max_dist, (lon_average, lat_average), max_coord

def producir_resultados(mmsi, list_coord):
    """
    Calcula e imprime para un barco en particular.

    :param mmsi: El barco.
    :param list_coord: La lista de coordenadas
    :return: Maxima distancia recorrida, coord promedio y coord más lejana
    """
    print("MMSI", mmsi, "(", len(list_coord), " puntos)")
    if len(list_coord) > 1:
        max_dist, average, fartest = max_distance(list_coord)
    else:
        max_dist = 0.0
        average = fartest = list_coord[0]

    print("    Maxima distancia =", max_dist)
    print()

    return max_dist, average, fartest


def calcula_distancia(csv_file):
    """
    Calcula la distancia para los barcos en un archivos.

    :param csv_file: El archivo con coord.
    """
    with open(csv_file) as csv_file:
        csv_reader = csv.reader(csv_file, delimiter=',')
        mmsi_act = None
        list_coord = []
        header = next(csv_reader, None)
        for mmsi, datetime, lon, lat in csv_reader:
            coord = (float(lon), float(lat))
            if mmsi_act is None:
                mmsi_act = mmsi
            if mmsi_act != mmsi:
                producir_resultados(mmsi_act, list_coord)
                mmsi_act = mmsi
                list_coord = []

            list_coord.append(coord)
        producir_resultados(mmsi_act, list_coord)

#
#   Revisar distancias usando
#   https://www.meridianoutpost.com/resources/etools/calculators/calculator-latitude-longitude-distance.php?
#
#   Ejemplo 1: Valores correctos son 49.79 miles / 80.13 kilometers / 43.24 nautical miles
#
print("Test para MMSI 209698000")
mmsi_209698000 = [(-1.89616501331329, 36.9659996032715),
                  (-0.959831655025482, 37.571834564209),
                  (-0.959831655025482, 37.571834564209),
                  (-0.959831655025482, 37.571834564209)]

max_dist, pos_average, pos_fartest = max_distance(mmsi_209698000)
print("    Distancia = ", max_dist, " millas naúticas.")

#
# Ejemplo 2: El barco no se mueve de su posición.
#
print("Test para MMSI 100000000")
mmsi_100000000 = [(-5.35431146621704, 36.1574516296387),
                  (-5.35431146621704, 36.1574516296387),
                  (-5.35431146621704, 36.1574516296387),
                  (-5.35431146621704, 36.1574516296387),
                  (-5.35431146621704, 36.1574516296387)]
max_dist, pos_average, pos_fartest = max_distance(mmsi_100000000)
print("    Distancia = ", max_dist, " millas naúticas.")

#
# Ejemplo 3: Pequeña distancia
# Valores correctos: 0.00 miles / 0.00 kilometers / 0.00 nautical miles
#
print("Test para MMSI 205834010")
mmsi_205834010 = [(-2.60475993156433 , 36.7570495605469),
                  (-2.60475826263428 , 36.7570571899414),
                  (-2.6047465801239 , 36.7570381164551),
                  (-2.60472822189331 , 36.7570648193359),
                  (-2.60472321510315 , 36.7570724487305),
                  (-2.60469508171082 , 36.7570648193359),
                  (-2.60473656654358 , 36.7570724487305),
                  (-2.60474991798401 , 36.7570610046387),
                  (-2.60470509529114 , 36.7570533752441)]
max_dist, pos_average, pos_fartest = max_distance(mmsi_205834010)
print("    Distancia = ", max_dist, " millas naúticas.")

print()
print("---- Proceso archivo ----")
calcula_distancia('barcos.csv')

Produce

Test para MMSI 209698000
    coord average = (-1.193914994597434, 37.420375823974624)
    coord fartest = (-1.89616501331329, 36.9659996032715)
    Distancia =  48.59234211469109  millas naúticas.
Test para MMSI 100000000
    coord average = (-5.35431146621704, 36.1574516296387)
    coord fartest = (-5.35431146621704, 36.1574516296387)
    Distancia =  0.0  millas naúticas.
Test para MMSI 205834010
    coord average = (-2.6047336525387244, 36.757059309217674)
    coord fartest = (-2.60469508171082, 36.7570648193359)
    Distancia =  0.0013938850605571295  millas naúticas.

---- Proceso archivo ----
MMSI      17620  ( 2418  puntos)
    coord average = (-6.180518489518966, 36.497049667699166)
    coord fartest = (-6.18065643310547, 36.4969139099121)
    Maxima distancia = 0.008632633913399173

MMSI      20233  ( 40  puntos)
    coord average = (-6.051580440998078, 35.45204601287842)
    coord fartest = (-6.05156850814819, 35.4520416259766)
    Maxima distancia = 0.0005171909496390382

MMSI    2424103  ( 1  puntos)
    Maxima distancia = 0.0

MMSI    2424105  ( 227  puntos)
    coord average = (-4.661880942693366, 35.21154432254741)
    coord fartest = (-4.66185998916626, 35.2115249633789)
    Maxima distancia = 0.0015238578744766443

...

MMSI  620069000  ( 8  puntos)
    coord average = (-0.30466312170028675, 35.85506343841553)
    coord fartest = (-0.304684996604919, 35.8550300598145)
    Maxima distancia = 0.002023965931953086

MMSI  636012425  ( 1  puntos)
    Maxima distancia = 0.0

MMSI  671222100  ( 62  puntos)
    coord average = (-2.466488011421696, 36.83163489064863)
    coord fartest = (-2.46646165847778, 36.8316688537598)
    Maxima distancia = 0.0022689686308574636
2
  • Gracias por tu respuesta! La idea sería sustituir las nubes de puntos que genera un barco "inmovil", es decir, que su posición apenas varíe en el tiempo, por un punto representativo. Quizás este punto pueda ser el centroide. Sobre los puntos "móviles", se considerarían de interés, por lo que a priori se dejan intactos. Un saludo! Commented el 3 mar. 2021 a las 22:05
  • Es decir, simplemente se busca eliminar posiciones de un mismo barco en la misma zona durante un tiempo, de manera que a la hora de usar métodos de clustering, no se genere un grupo con ese barco solo Commented el 3 mar. 2021 a las 22:07
0

Una idea inspirada en el comentario de CandidMoe.

Cargamos el csv en un dataframe:

import pandas as pd

df = pd.read_csv("trozo2.csv", sep="\s*,\s*")

Si ahora agrupamos todas las filas que corresponden a un mismo barco (agrupar por "mmsi") y a cada grupo le calculamos la desviación estándar de las columnas "latitud" y "longitud", los que tengan mayor desviación estándar serán los que más se hayan movido (y a la inversa los de menor, los que menos).

Esto podríamos hacerlo con:

df.groupby("mmsi").std()

Sin embargo sería más útil tener un solo número para la desviación estándar, en vez de dos (uno para latitud y otro para longitud), por lo que se me ocurre considerar cada desviación como un vector con dos componentes (desviación de la latitud y de la longitud), y calcularíamos entonces la longitud del vector usando el teorema de Pitágoras (en esencia: raiz cuadrada de la suma de los cuadrados de ambas desviaciones).

Esto puede hacerse de forma relativamente simple con math.hypot():

import math
std = df.groupby("mmsi").std().apply(
        lambda x: math.hypot(x.longitude, x.latitude), axis=1).fillna(0)

El fillna(0) es para llenar con 0 los casos que hayan dado NaN. Estos serán los casos para los cuales el barco ha aparecido una sola vez (en ese caso no se puede calcular la desviación típica). Les asigno 0 porque en ese caso los considero inmóviles. Si quisieras considerarles móviles podrías asignarles un valor arbitrario, por ejemplo 2.

Como consecuencia en std tenemos ahora algo así:

mmsi
17620        0.000036
20233        0.000005
2424103      0.000000
2424105      0.000012
2424106      0.000018
               ...   
605066435    0.000034
605116715    0.000022
620069000    0.000022
636012425    0.000000
671222100    0.000025

Añadimos estos datos al dataframe original:

df = df.set_index("mmsi").sort_index()  # Para que df tenga el mismo indice que std
df["std"] = std

El dataframe tiene ahora una columna std que de algún modo nos dice cuánto se ha movido cada barco (y tiene el mismo valor repetido para todas las filas en que aparece ese barco). Por ejemplo, df empieza ahora así:

>>> df.head()
                  msgtimestamp  longitude   latitude       std
mmsi                                                          
17620  2019-12-11 05:54:22.954  -6.180553  36.497028  0.000036
17620  2019-12-11 05:56:03.025  -6.180430  36.497025  0.000036
17620  2019-12-11 05:57:15.299  -6.180459  36.497036  0.000036
17620  2019-12-11 05:58:23.291  -6.180483  36.497021  0.000036
17620  2019-12-11 06:01:44.753  -6.180521  36.497051  0.000036

Estos números en la columna std no nos dicen mucho, pero si pedimos la descripción estadística de estos datos con std.describe() sale:

>>> std.describe()
count    222.000000
mean       0.017102
std        0.093680
min        0.000000
25%        0.000000
50%        0.000025
75%        0.000569
max        1.116976

Aquí hay información relevante. 222 es el número de barcos diferentes. Entre ellos los hay que no se han movido nada, pues el mínimo sale 0.0000 Y para el que más se ha movido la desviación ha salido 1.116976. En los percentiles vemos que el 75% de los barcos ha dado un valor de std inferior a 0.000569

Imagino que estos datos podrían tener más sentido para ti que para mí. El caso es que a partir de ellos habría que elegir un cierto umbral para std por debajo del cual consideraremos que "el barco no se ha movido" y por encima del cual "sí se ha movido".

Yo tomaré un valor relativamente arbitrario, como 0.001. Entonces:

df["inmovil"] = df["std"] < 0.001

Si tomamos una muestra de 10 filas al azar del dataframe resultante podemos ver que aparecen barcos marcados como inmóviles, y otros que no:

>>> df.sample(10)
                      msgtimestamp  longitude   latitude       std  inmovil
mmsi                                                                       
17620       2019-12-11 06:31:28.07  -6.180553  36.497002  0.000036     True
224225750  2019-12-11 23:04:35.397  -2.606215  36.754967  0.000036     True
17620      2019-12-12 03:38:34.743  -6.180517  36.497028  0.000036     True
209405000   2019-12-12 13:48:26.56  -5.517715  35.878517  0.132696    False
2424140    2019-12-11 15:15:05.224  -2.921977  35.270615  0.000012     True
224363830  2019-12-11 12:43:09.895  -5.933078  36.186394  0.000024     True
242060100  2019-12-11 09:39:08.106  -5.808740  35.792912  0.000210     True
224363830  2019-12-12 15:30:49.966  -5.933052  36.186367  0.000024     True
224062640  2019-12-11 05:47:20.182  -2.482865  36.835087  0.003583    False
224363830  2019-12-12 01:54:23.885  -5.933045  36.186394  0.000024     True

A partir de aqui ya puedes usar esa columna para filtrar. Por ejemplo vemos que de las 24232 filas (len(df)) de tu tabla, la mayoría correspondería a barcos "inmóviles" según el criterio usado:

>>> len(df[df["inmovil"]])
19356
1
  • Muchas gracias! Aunque no es la solución exacta a lo que preguntaba, me aporta información de gran utilidad. SI me lo permites, enlazo otra pregunta... A partir de este punto, es posible sustituir los barcos "inmóviles" por algún punto representativo?? (por ejemplo, su centroide). Insisto en que soy bastante novato en esto, por lo que reitero mis disculpas si pregunto algo básico... Gracias!! :) Commented el 3 mar. 2021 a las 21:58

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.