Theil-Sen-Schätzerfunktion in T-SQL

Antworten:

9

Ich habe gelogen, als ich sagte, ich kann es nicht in SQL umcodieren. Ich war einfach zu faul. Hier ist der Code mit einem Verwendungsbeispiel.

Der Code basiert auf einer TheiSen- Perl-Bibliothek unter Verwendung von QuickMedian . Definieren wir einen neuen Tabellentyp, um unsere Daten einfach an die Prozedur zu übergeben.

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

Bitte beachten Sie die ID-Spalte, die hier wichtig ist, da unsere Lösung die Anweisung CROSS APPLY verwendet, um eine korrekte Interpretation der inneren Schleife in TheilSen.pm zu erreichen.

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];

Wir benötigen außerdem einen neuen Datentyp zum Speichern eines Arrays realer Typwerte.

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

Hier ist das Funktion f_QuickMedian , die den Median für ein bestimmtes Array zurückgibt . Der Kredit dafür geht an 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

Und die p_TheilSen- Schätzer:

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

Beispiel:

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

Kehrt zurück:

m = 52.7079
c = -448.4853

Nur zum Vergleich gibt die Perl-Version die folgenden Werte für denselben Datensatz zurück:

m = 52.7078651685394
c = -448.484943820225

Ich verwende den TheilSen-Schätzer, um die DaysToFill-Metrik für Dateisysteme zu berechnen. Genießen!

kafe
quelle
1

Ich habe auch nach T-SQL, Oracle und Servern im Allgemeinen gesucht (zu komplex, um in reinem SQL geschrieben zu werden).

In Allerdings könnten Sie interessieren diese (ein wissenschaftliches / Statistik - Paket für Python). Der Algorithmus ist dort auch und in Python implementiert. Python ist eine Sprache, die Menschen im Gegensatz zu Perl zumindest leicht verstehen können.

Ihre Frage hat mich fasziniert und ich habe mich umgegraben. Es gibt C- und C ++ - Bibliotheken, die diesen Algorithmus enthalten - und er ist auch in einigen R-Paketen verfügbar. Und auch @srutzkys Beitrag sieht interessant aus.

+1 für eine interessante Frage übrigens - und willkommen im Forum :-)

Vérace
quelle