(Python) Estimating Regressionsparameter Konfidenzintervall mit scikits Bootstrap

Question

Ich habe gerade begonnen, durch scikits verfügbar ein schönes Bootstrapping Paket ausprobieren: https://github.com/cgevans/scikits-bootstrap

aber ich habe ein Problem festgestellt, wenn sie versuchen Konfidenzintervall abzuschätzen für der Korrelationskoeffizient aus der linearen Regression. Die zurückgegebenen Konfidenzintervalle liegen vollständig außerhalb des Bereichs der ursprünglichen Statistik. Hier

ist der Code:

import numpy as np 
from scipy import stats 
import bootstrap as boot 

np.random.seed(0) 
x  = np.arange(10) 
y  = 10 + 1.5*x + 2*np.random.randn(10) 
r0 = stats.linregress(x, y)[2] 

def my_function(y): 
    return stats.linregress(x, y)[2] 

ci = boot.ci(y, statfunction=my_function, alpha=0.05, n_samples=1000, method='pi') 

dies ein Ergebnis von ci ergibt = [-0,605, 0,644], aber die ursprüngliche Statistik ist r0 = 0,894.

Ich habe das in R versucht und es scheint dort gut zu funktionieren: Die Ci überspannt R0 wie erwartet.

Bitte helfen!

Answer 1

Können Sie Ihren R-Code angeben? Ich würde interessiert sein zu wissen, wie dies in R. behandelt wird.

Das Problem hier ist, dass Sie nur y zu boot.ci übergeben, aber jedes Mal, wenn es my_function ausgeführt wird, verwendet es die gesamte, ursprüngliche x (Beachten Sie das Fehlen von x-Eingabe für my_function). Bootstrapping wendet die Statistikfunktion auf neu abgetastete Daten an. Wenn Sie also Ihre Statistikfunktion mit dem Original x und einem Beispiel y anwenden, haben Sie ein unsinniges Ergebnis. Aus diesem Grund funktioniert die BCA-Methode überhaupt nicht: Sie kann Ihre Statistikfunktion nicht auf Jackknife-Samples anwenden, die nicht die gleiche Anzahl an Elementen haben. Wenn Sie mehrere 1D-Arrays an Ihre Statistikfunktion übergeben möchten, können Sie mehrere Spalten verwenden: xy = vstack((x,y)).T würde funktionieren und dann eine Statistikfunktion verwenden, die Daten aus diesen Spalten übernimmt:

def my_function(xysample): 
    return stats.linregress(xysample[:,0], xysample[:,1])[2] 

Alternativ, wenn Sie mit Ihren Daten überhaupt in Unordnung vermeiden wollten, könnten Sie eine Funktion definieren, die auf Indizes arbeitet, und dann nur Indizes passieren boot.ci:

def my_function2(i): 
    return stats.linregress(x[i], y[i])[2] 

boot.ci(np.arange(len(x)), statfunction=my_function2, alpha=0.05, n_samples=1000, method='pi') 

Beachten Sie, dass in jedem dieser Fälle BCA funktioniert, Sie können also auch method = 'bca' verwenden, es sei denn, Sie möchten wirklich prozentuale Intervalle verwenden; BCA ist ziemlich immer besser.

Ich bin mir bewusst, dass diese beiden Methoden nicht ideal sind. Ehrlich, ich hatte noch nie eine Notwendigkeit, mehrere Arrays wie diese an meine Statfunktion zu übergeben, und die meisten Leute benutzen wahrscheinlich mean als ihre Statfunktion. Ich denke, die beste Idee hier könnte sein, Listen von gleichgroßen [0] Arrays zu übergeben, zB boot.ci([x,y],...), und dann alle diese gleichzeitig abzutasten und sie alle als separate Argumente an die statfunction zu übergeben. In diesem Fall könnten Sie einfach eine my_function(x,y) haben. Ich werde sehen, ob ich das machen kann, aber wenn du mir deinen R-Code zeigen kannst, wäre das großartig, denn ich würde gerne sehen, ob es einen besseren Weg gibt, damit umzugehen.

---

Update:

In der neuesten Version von scikits.bootstrap (v0.3.1), ein Tupel von Arrays zur Verfügung gestellt werden kann, und Proben von ihnen werden als separate Argumente statfunction weitergegeben werden. Darüber hinaus kann statfunction Array-Ausgaben bereitstellen, und Konfidenzintervalle werden für jeden Punkt in der Ausgabe berechnet.Somit ist dies jetzt sehr einfach zu tun. Im Folgenden wird geben Konfidenzintervalle für jede Ausgabe von linregress:

cis = boot.ci((x,y), statfunction=stats.linregress) 

cis[:,2] in diesem Fall wird das gewünschte Vertrauensintervall sein.