Función de estimador Theil-Sen en T-SQL

Respuestas:

9

Estaba mintiendo cuando dije, no puedo volver a codificarlo en SQL. Era demasiado vago. Aquí está el código con un ejemplo de uso.

El Código se basa en una biblioteca perl TheiSen , usando QuickMedian . Definamos un nuevo tipo de tabla para pasar fácilmente nuestros datos al procedimiento.

CREATE TYPE dbo.TheilSenInputDataTableType AS TABLE 
(
    ID INT IDENTITY(1,1),
    x REAL, 
    y REAL
)

Tenga en cuenta la columna ID, que es importante aquí, ya que nuestra solución utiliza la declaración CROSS APPLY para lograr una interpretación correcta del bucle interno que se encuentra en TheilSen.pm.

my ($x1,$x2,$y1,$y2);
foreach my $i(0 .. $n-2){
    $y1 = $y->[$i];
    $x1 = $x->[$i];
    foreach my $j($i+1 .. $n-1){
        $y2 = $y->[$j];
        $x2 = $x->[$j];

También necesitaremos un nuevo tipo de datos para almacenar una matriz de valores de tipo real.

CREATE TYPE [dbo].[RealArray] AS TABLE(
    [val] [real] NULL
)

Aquí está la función f_QuickMedian , que devuelve la mediana para una matriz dada. El crédito para este va a Itzik Ben-Gan .

CREATE FUNCTION [dbo].[f_QuickMedian](@RealArray RealArray READONLY)
RETURNS REAL
AS
BEGIN
    DECLARE @Median REAL;
    DECLARE @QMedian REAL;

    SELECT @Median = AVG(1.0 * val)
    FROM
    (
        SELECT o.val, rn = ROW_NUMBER() OVER (ORDER BY o.val), c.c
        FROM @RealArray AS o
        CROSS JOIN (SELECT c = COUNT(*) FROM @RealArray) AS c
    ) AS x
    WHERE rn IN ((c + 1)/2, (c + 2)/2);

    SELECT TOP 1 @QMedian = val FROM @RealArray
    ORDER BY ABS(val - @Median) ASC, val DESC

    RETURN @QMedian
END

Y el estimador p_TheilSen :

CREATE PROCEDURE [dbo].[p_TheilSen](
      @TheilSenInput TheilSenInputDataTableType READONLY
    , @m Real OUTPUT
    , @c Real OUTPUT
)
AS
BEGIN
    DECLARE 
        @m_arr RealArray
      , @c_arr RealArray;       

    INSERT INTO @m_arr
        SELECT m
        FROM 
        (
            SELECT  
                t1.x as x1
                , t1.y as y1
                , t2o.x as x2
                , t2o.y as y2
                , t2o.y-t1.y as [y2 - y1]
                , t2o.x-t1.x as [x2 - x1]
                , CASE WHEN (t2o.x <> t1.x) THEN  CAST((t2o.y-t1.y) AS Real)/(t2o.x-t1.x) ELSE NULL END AS [($y2-$y1)/($x2-$x1)]
                , CASE WHEN t1.y = t2o.y THEN 0
                  ELSE
                    CASE WHEN t1.x = t2o.x THEN NULL
                        ELSE 
                        -- push @M, ($y2-$y1)/($x2-$x1);
                        CAST((t2o.y-t1.y) AS Real)/(t2o.x-t1.x)
                    END
                  END as m
            FROM @TheilSenInput t1
            CROSS APPLY
                    (
                    SELECT  t2.x, t2.y
                    FROM    @TheilSenInput t2
                    WHERE   t2.ID > t1.ID
                     ) t2o
        ) t
        WHERE m IS NOT NULL 

    SELECT @m = dbo.f_QuickMedian(@m_arr)

    INSERT INTO @c_arr
        SELECT y - (@m * x)
            FROM @TheilSenInput

    SELECT @c = dbo.f_QuickMedian(@c_arr)

END

Ejemplo:

DECLARE 
      @in TheilSenInputDataTableType
    , @m Real
    , @c Real

INSERT INTO @in(x,y) VALUES (10.79,118.99)
INSERT INTO @in(x,y) VALUES (10.8,120.76)
INSERT INTO @in(x,y) VALUES (10.86,122.71)
INSERT INTO @in(x,y) VALUES (10.93,125.48)
INSERT INTO @in(x,y) VALUES (10.99,127.31)
INSERT INTO @in(x,y) VALUES (10.96,130.06)
INSERT INTO @in(x,y) VALUES (10.98,132.41)
INSERT INTO @in(x,y) VALUES (11.03,135.89)
INSERT INTO @in(x,y) VALUES (11.08,139.02)
INSERT INTO @in(x,y) VALUES (11.1,140.25)
INSERT INTO @in(x,y) VALUES (11.19,145.61)
INSERT INTO @in(x,y) VALUES (11.25,153.45)
INSERT INTO @in(x,y) VALUES (11.4,158.03)
INSERT INTO @in(x,y) VALUES (11.61,162.72)
INSERT INTO @in(x,y) VALUES (11.69,167.67)
INSERT INTO @in(x,y) VALUES (11.91,172.86)
INSERT INTO @in(x,y) VALUES (12.07,177.52)
INSERT INTO @in(x,y) VALUES (12.32,182.09)


EXEC p_TheilSen @in, @m = @m OUTPUT, @c = @c OUTPUT

SELECT @m
SELECT @c

Devoluciones:

m = 52.7079
c = -448.4853

Solo para una comparación, la versión perl devuelve los siguientes valores para el mismo conjunto de datos:

m = 52.7078651685394
c = -448.484943820225

Utilizo el estimador TheilSen para calcular la métrica DaysToFill para sistemas de archivos. ¡Disfrutar!

kafe
fuente
1

También verifiqué, para T-SQL, Oracle y servidores en general (demasiado complejo para ser escrito en SQL puro).

Sin embargo, puede interesarle esto (un paquete científico / estadístico para Python). El algoritmo se implementa allí también y en Python. Python es un lenguaje que los humanos tienen al menos alguna posibilidad de poder entender, a diferencia de Perl.

Su pregunta me intrigó y busqué. Hay bibliotecas C y C ++ que contienen este algoritmo, y también está disponible en algunos paquetes R. Y la publicación de @srutzky también parece interesante.

+1 para una pregunta interesante BTW - y bienvenido al foro :-)

Vérace
fuente