STATA-FORUM.DE

[dpark]

Hallo zusammen,

ich habe folgendes Problem. Ich habe drei Schleifen, die jeweils 18 Regressionen generieren. In der ersten Schleife gibt es nur die dependent variable und eine independent variable. Die zweite fügt einen Interaction term mit einer Dummy hinzu und die dritte Schleife nochmal additiv die Dummy Variable.

Da das zugrunde liegende Model ein probit model ist und ich robuste SE benötige, verwende ich die funktion nwest probit depvar indepvar, lag(x). Mein Code momentan sieht wie folgt aus:

Code: Alles auswählen

forvalues i=1/8 {
nwest probit depvar var1_`i', lag(`i')
nwest probit depvar var2_`i', lag(`i')
nwest probit depvar var3_`i', lag(`i')
nwest probit depvar var4_`i', lag(`i')
nwest probit depvar var5_`i', lag(`i')
nwest probit depvar var6_`i', lag(`i')
}
forvalues j=1/8 {
nwest probit depvar var1_`j' var1_wd`j', lag(`j')
nwest probit depvar var2_`j' var2_wd`j', lag(`j')
nwest probit depvar var3_`j' var3_wd`j', lag(`j')
nwest probit depvar var4_`j' var4_wd`j', lag(`j')
nwest probit depvar var5_`j' var5_wd`j', lag(`j')
nwest probit depvar var6_`j' var6_wd`j', lag(`j')
}
forvalues k=1/8 {
nwest probit depvar var1_`k' var1_wd`k' dummy`k', lag(`k')
nwest probit depvar var2_`k' var2_wd`k' dummy`k', lag(`k')
nwest probit depvar var3_`k' var3_wd`k' dummy`k', lag(`k')
nwest probit depvar var4_`k' var4_wd`k' dummy`k', lag(`k')
nwest probit depvar var5_`k' var5_wd`k' dummy`k', lag(`k')
nwest probit depvar var6_`k' var6_wd`k' dummy`k', lag(`k')
}


Die Variable varX_wd ist der Interaction term, den ich manuell ausrechnen musste, da nwest keine operators erlaubt.

Nun will ich die berechneten t-stats für die jeweiligen Variablen gerne extrahieren, d.h. nach der ersten Schleife einen Befehl einfügen, der mir die t-Stat Werte in eine Zeile einer Matrix speichert.
Für die zweite Spalte entsprechend eine 2x8 Matrix, wobei die t-Stats der independentvar in die erste Zeile sollen und die t-Stats des Interaction terms in die zweite Zeile.
Analog alles für die dritte Schleife.

Kann mir einer einen Tipp geben, wie ich dies am einfachsten löse? Ich hab es bereits einmal für eine andere depVar gemacht und immer händisch die Werte rauskopiert. Das ist aber sehr mühselig. Hatte vorher mit Eviews gearbeitet, was in der Programmierung deutlich einfacher ist. Leider kann Eviews keine robusten SE nach Newey-West für probit berechnen.

Vielen Dank im Voraus!

[daniel]

Den Inhalt zu kommentieren, fehlt mir die Zeit.

Der t-Wert berechnet sich bekanntlich als



Nach der Regression sind die Koeffizienten in der Matrix e(b) zu finden. Die Covarianzmatrix ist in e(V) hinterlegt.

Code: Alles auswählen

sysuse auto ,clear
reg price mpg foreign
mat b = e(b)
mat V = e(V)
di "t(mpg) = " b[1, 1]/sqrt(V[1, 1])


Falls diese Antwort zufriedenstellend ist, bitte daran denken den thread auf Statalist (http://www.stata.com/statalist/archive/ ... 00798.html) mit einer cross-reference hierher zu schließen, wie kürzlich von Nick Cox vorgeschlagen (http://www.stata.com/statalist/archive/ ... 00817.html).

btw. Du solltest auf Statalist erklären, wo Du user-written commands her hast. Das ist auch hier eine gute Idee.

STB-39 sg72 Newey-West standard errors for probit, logit, and poisson models
(help nwest if installed) . . . . . . . . . . . . . . . . . J. Hardin
9/97 pp32--35; STB Reprints Vol 7, pp.182--186

[dpark]

Ich sitze immer noch dran. Sobald ich das Problem aus der Welt geschafft hab, mach ich auf Statalist mit cross reference zu.

Bin noch am basteln was das Zusammenfassen der Ergebnisse in einer übersichtlichen Matrix angeht.

EDIT: Also hier ein Ausschnitt meines Codes zur Berechnung der Pseudo-R² (gemäß Estrella-Mishkin 1998, eine leicht abgewandelte Form des McFaddden R²)

Code: Alles auswählen

*Calculation of Pseudo-R² for Variables alone
forvalues a=1/8 {
probit ptinrec spread`a'
mat m_spread`a'=1-(e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
probit ptinrec dlog_rm0_`a'
mat m_rm0_`a'=1-(e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
probit ptinrec dlog_rgdp`a'
mat m_rgdp`a'=1-(e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
probit ptinrec dlog_cpi`a'
mat m_cpi`a'=1-(e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
probit ptinrec dlog_spcomp`a'
mat m_spcomp`a'=1-(e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
}
*Values of Pseudo-R² into Matrices
mat m_spread_full=(m_spread1,m_spread2,m_spread3,m_spread4,m_spread5,m_spread6,m_spread7,m_spread8)
mat m_rm0_full=(m_rm0_1,m_rm0_2,m_rm0_3,m_rm0_4,m_rm0_5,m_rm0_6,m_rm0_7,m_rm0_8)
mat m_rgdp_full=(m_rgdp1,m_rgdp2,m_rgdp3,m_rgdp4,m_rgdp5,m_rgdp6,m_rgdp7,m_rgdp8)
mat m_cpi_full=(m_cpi1,m_cpi2,m_cpi3,m_cpi4,m_cpi5,m_cpi6,m_cpi7,m_cpi8)
mat m_spcomp_full=(m_spcomp1,m_spcomp2,m_spcomp3,m_spcomp4,m_spcomp5,m_spcomp6,m_spcomp7,m_spcomp8)
*Final Matrix for all Variables alone
mat single=(m_spread_full\m_rm0_full\m_rgdp_full\m_cpi_full\m_spcomp_full)


Kann ich das noch irgendwie vereinfacht programmieren? Ich hab mit Stata erst am Sonntag angefangen, daher bitte ich um Entschuldigung für diese eventuell Einsteigerfrage.

[daniel]

Hier Dein code, wie ich ihn schreiben würde (nicht ganz -- ich würde das noch allgemeiner schreiben)

Code: Alles auswählen

// get rid of underscore in m_rm0_#
ren (m_rm0_#)(m_rm0#)

// get all variable stubs into local
loc stubs m_spread m_rm0 m_rgdp m_cpi m_spcomp

// get variable stubs into local 1, 2, ...
token `stub'

// set up results matrix
mat M = J(5, 8, .)
mat rownam M = `stubs'

// run probit models and fill in r2
forv r = 1/5 {
   forv c = 1/8 {
      qui prob ptinrec ``r''`c'
      mat M[`r', `c'] = 1 - (e(ll)/e(ll_0))^(-(2/e(N))*e(ll_0))
   }
}


[dpark]

Hallo Daniel.

Vielen Dank für Deine Hilfe. Auch wenn ich zugegebenermaßen Deinen Code nicht zu 100 Prozent verstanden habe (konzeptionell schon, nur das mit tokenize noch nicht ganz), ist mein Code jetzt fertig. Ist immer noch kompliziert und nicht so schön kompakt wie Deiner, aber er kann zumindest alles ausrechnen, was ich brauche.

Eine letzte Frage hätte ich noch. Ist nur indirekt im Bezug auf Syntax, aber ich wollte nicht einen extra Thread aufmachen.

Ich mache eine out-of-sample prediction. Als Ergebnis will ich zwei Dinge haben: Graph mit der predicted probability und den Wert für das Pseudo-R², wie schon oben eingebaut.

Den Graphen krieg ich ja sehr leicht mit

Code: Alles auswählen

local start 1993q2
probit depvar indepvar if date<tq(`start')
predict depvar_forecasted if date>=tq(`start'), p
twoway (line depvar_forecasted date) if date>tq(`start'), ylabel(0(0.2)1)


Allerdings macht mir das Pseudo-R² Sorgen als das ich nicht ganz damit klarkomme. Der Wert berechnet sich (ich ignoriere den Exponenten erstmal) als Differenz zwischen der unrestricted Log Likelihood und der restricted Log Likelihood geteilt durch die restricted Log Likelihood.

Den Wert für die unrestricted Log Likelihood kann ich ja einfach mit den predicted probabilities als Output des "predict" Befehls ausrechnen. Wie aber komme ich an die Werte für das restricted Log Likelihood ran, also die predicted probabilities für ein constant-only model? Eine solche Einschränkung habe ich für den "predict" Befehl leider nicht gefunden.

EDIT:
Okay, ich hab jetzt verstanden, wie ich manuell rankommen könnte. Werde es dann mal versuchen.

[daniel]

nur das mit tokenize noch nicht ganz

Steht ja alles in den help files etc., aber zur kurzen Erläuterung. -tokenize- spaltet eine Zeichenkette in ihre Bestandteile (tokens), wobei die einzelenen Bestandteile (per default) durch Leerzeiche voneinander getrennt sind. Das Ergebnis wird in den local macros 1, 2, ... gespeichert. Der "Trick" an der Sache ist, dass local genestet sein dürfen (da gibt es Grenzen, weiß nicht mehr wo genau die liegen). Hier nochmal ein Minibsp, dass die Idee (etwas auf die Spitze getrieben) verdeutlicht

Code: Alles auswählen

loc foobar first second third
loc first foo
loc second bar
loc third foobar

token `foobar'
forv j = 1/3 {
   di
   di _col(2) as txt "{char 96}j{char 39} is now " as res "`j' " ///
   as txt "and {char 96}{char 96}j{char 39}{char 39} is " ///
   "{char 96}`j'{char 39} and therefore "  as res "``j''"
   if (`j' == 3) {
      di _col(2) as txt "and since " ///
      "{char 96}{char 96}{char 96}j{char 39}{char 39}{char 39} " ///
      "is {char 96}{char 96}`j'{char 39}{char 39} " ///
      "which expands to {char 96}``j''{char 39} " ///
      "we are back at " as res "```j'''"
   }
}


Wie aber komme ich an die Werte für das restricted Log Likelihood ran, also die predicted probabilities für ein constant-only model?

Am leichtesten, indem Du ein constant-only model schätzt, und im Anschluss -predict-est. Möglich, dass Du mit -margins- auch weiterkommst, ist aber komplexer für den Einstieg.

[dpark]

Super, vielen Dank für das Beispiel. Sobald ich mit meinem Code durch bin, schau ich mir das mal genau an. An sich finde ich Stata super, nur für Einsteiger ist es nicht zwingend die beste Wahl

Ich hatte gedacht, mit dem constant-only ginge es eigentlich relativ einfach. Aber zB.

Ich regressiere probit depvar indepvar. Im Output meldet er ja die Iterations für das LL. Ich habe es so verstanden (und im Internet gelesen), dass die 0-te Iteration gleich dem restricted LL ist. Leider kann ich das nicht abgreifen, denn dann wäre es extrem einfach.

Wenn ich aber probit depvar einfach so regressiere, dann komme ich auf ein ganz anderes LL als die 0-te Iteration. Ist meine Annahme falsch, dass die 0-te Iteration gleich dem restricted LL ist (wiewohl sich dieser Wert mit dem Wert der von Eviews für selbige Regression deckt) oder wo könnte der Fehler liegen?

Auf jeden Fall nochmal vielen Dank für die tolle Hilfe!

[daniel]

An sich finde ich Stata super, nur für Einsteiger ist es nicht zwingend die beste Wahl

Ich kenne sonst nur SPSS und R und muss sagen, da ziehe ich Stata, gerade als Anfänger, beiden deutlich vor. Was ich bei Stata so einfach finde, ist die durchgängig konstante (einfache) Syntax. Als Sprache betrachtet ist also die Grammatik exterm intuitiv und einfach zu erlernen, und man muss sich eben im Laufe der Zeit nur noch die Vokabeln draufschaffen.

Leider kann ich das nicht abgreifen, denn dann wäre es extrem einfach.

Doch das geht, sorry, dass ich daran nicht früher gedacht habe -- allerdings hättest Du das auch im help file nachlesen können. Die log-likelihood des constant-only models ist in e(ll_0) hinterlegt. Das default Pseudo-R Square (ich meine, Stata verwendet McFadden) ist direkt in e(r2_p) hinterlegt.

[dpark]

Hallo Daniel,

das mit dem Abreifen via e(ll_0) kenne ich ja. Allerdings geht das ja nicht, da ich für meine prediction keine normale Regression laufen lasse. Mein Ziel sind rolling-sample out-of-sample forecasts. Dazu berechne ich die predicted probabilities meiner depvar sukzessive mit den neuen Koeffizienten der "sich rollenden" probit-regression, indem ich den wert der standard normal cumulative distribution meiner entsprechenden linear-kombination mit der "normal(z)" funktion ausgeben lasse. Meine unrestricted LL erhalte ich dann durch das Produkt der Wahrscheinlichkeiten (bzw. der Gegenwahrscheinlihckeit für den Fall dass meine DepVar 0 ist.).

An dieser Stelle wollte ich fragen, ob es eine Funktion gibt, das Produkt der Elemente einer Spalte in einer Matrix zu berechnen. Ich hab dazu einfach nichts gefunden. Ich wollte es mit einem Trick versuchen und das Produkt in meine Schleife einbauen, sodass bei der Berechnung jeder neuen Wahrscheinlichkeit diese ausgerechnete WK mit dem Produkt der bisherigen WKs multipliziert wird. Aber irgendwie kriege ich das einfach nicht gebacken.

Code: Alles auswählen

tsset date

local horizon "1qtr 2qtrs 3qtrs 4qtrs 5qtrs 6qtrs 7qtrs 8qtrs"
local indepvarlist "spread rm0 rgdp cpi spcomp"

*Out-of-sample forecast
local null 236
local start 360
local max=453-`start'-1

foreach indepvar in `indepvarlist' {
mat prob01_ptinrec_`indepvar'=J(`max', 8, .)

forvalues i=1/8 {
forvalues c=1/`max' {
local current=`start'+`c'
probit ptinrec `indepvar'`i' in `null'/`current'
mat a=e(b)
if ptinrec[`current'+1]>0 {
mat prob01_ptinrec_`indepvar'[`c',`i'] = normal(-(a[1,2]+a[1,1]*`indepvar'`i'[`current']))
}
else {
mat prob01_ptinrec_`indepvar'[`c',`i'] = 1 - normal(-(a[1,2]+a[1,1]*`indepvar'`i'[`current']))
}
}
}
mat colnames prob01_ptinrec_`indepvar' = `horizon'
}

Da ist das Produkt noch nicht drin, weil ich es einfach nicht hingekriegt hab. Notfalls kann ich aber einfach die WKs in eine Excel exportieren und die unrestricted Log Likelihood extern berechnen.

Jetzt bleibt nur noch das Problem der restricted LL. Denn dazu müsste ich ja den Koeffizienten bzw die Konstante des restricted Models haben und damit die entsprechenden WKs ausrechnen. An sich ist das kein Problem. Allerdings komme ich immer auf ganz andere Werte, wenn ich es testweise für meine Estimation period ausrechne, d.h. ich lasse eine normale probit regression laufen und extrahiere die restricted LL. Dann lasse ich eine constant-only regression laufen und schau mir die LL an, die aber ganz anders ist als die restricted LL von vorher. Hab ich hier einfach einen Denkfehler?

PS:
Wenns die Funktion gäbe, würde ich dir jetzt mindestens 3 Tage lang nen Kaffee ausgeben

[daniel]

Mein Ziel sind rolling-sample out-of-sample forecasts.

Schau Dir mal -help rolling- an. Erleichtert evtl. die Sache.

An dieser Stelle wollte ich fragen, ob es eine Funktion gibt, das Produkt der Elemente einer Spalte in einer Matrix zu berechnen.

Ist mir nicht bekannt, aber das lässt sich leicht (in Mata) programmieren.

Code: Alles auswählen

vers 12.1
m :
real matrix mfcolprdc(real matrix x)
{
   real matrix p
   
   p = x[1, .]
   for (i = 2; i <= rows(x); ++i) {
      p = p :* x[i, .]
   }
   return(p)
}
end

// check how that works
mat foo = (1\ 2\ 3), (4\ 5\ 6)
matlist foo

// calcualte column products
m : mfcolprdc(st_matrix("foo"))

// Put the result into a new Stata matrix bar
m : st_matrix("bar", mfcolprdc(st_matrix("foo")))

matlist bar