[Rsiena-commits] r89 - in pkg/RSienaTest: R doc man

Sat May 29 23:10:40 CEST 2010

Author: jalospinoso
Date: 2010-05-29 23:10:40 +0200 (Sat, 29 May 2010)
New Revision: 89

Modified:
   pkg/RSienaTest/R/sienaTimeTest.r
   pkg/RSienaTest/R/sienaeffects.r
   pkg/RSienaTest/doc/RSiena.bib
   pkg/RSienaTest/doc/s_man400.tex
   pkg/RSienaTest/man/includeEffects.Rd
   pkg/RSienaTest/man/sienaTimeTest.Rd
Log:
Added functionality to sienaTimeTest, updated the manual.

Modified: pkg/RSienaTest/R/sienaTimeTest.r
===================================================================

--- pkg/RSienaTest/R/sienaTimeTest.r	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/R/sienaTimeTest.r	2010-05-29 21:10:40 UTC (rev 89)
@@ -10,7 +10,7 @@
 ## * sienaTimeFix, which is called to set up time dummy interacted effects.
 ## ****************************************************************************/
 ##@sienaTimeTest siena07 Does test for time homogeneity of effects
-sienaTimeTest <- function (sienaFit, effects=NULL)
+sienaTimeTest <- function (sienaFit, effects=NULL, condition=FALSE)
 {
 	observations <- sienaFit$f$observations
 	# There must be more than 2 observations to do a time test!
@@ -178,6 +178,7 @@
 				## these operations dont relate directly to the scores
 				dummyByEffect[i, j]=inc
 				dummyProps$shortName[inc] <- sienaFit$effects[-escreen,]$shortName[i]
+				dummyProps$interaction1[inc] <- sienaFit$effects[-escreen,]$interaction1[i]
 				dummyProps$type[inc] <- sienaFit$effects[-escreen,]$type[i]
 				dummyProps$period[inc] <- j + 1
 			}
@@ -192,17 +193,28 @@
 	doTests <- c(rep(FALSE, nEffects), rep(TRUE, nDummies))
 	jointTest <- ScoreTest(nTotalEffects, D, sigma, fra, doTests, maxlike=FALSE)
 	jointTestP <- 1 - pchisq(jointTest$testresOverall, nDummies)
-	individualTest <- jointTest$testresulto[1:nDummies]
+	if (! condition) {
+		individualTest <- jointTest$testresulto[1:nDummies]
+	} else {
+		individualTest <- sapply(1:nDummies, function (i)
+			{ doTests <- rep(FALSE, nEffects + nDummies)
+				doTests[nDummies+i] <- TRUE
+				test <- ScoreTest(nTotalEffects, D, sigma, fra, doTests, FALSE)
+				test$testresulto[1]
+			})
+	}
 	individualTestP <- 2 * (1-pnorm(abs(individualTest))[1:nDummies])
 	rownames(jointTestP) <- c("Joint Significant Test")
 	colnames(jointTestP) <- c("p-Val")
 	thetaOneStep <- c(sienaFit$theta[-c(dscreen,rscreen,escreen)], rep(0, nDummies)) +
-	jointTest$oneStep
+			jointTest$oneStep
 	effectTest <- sapply(1:length(indBaseEffects), function (i)
 					{
 						 doTests <- rep(FALSE, nEffects + nDummies)
 						 tmp <- which(dummyProps$shortName ==
-									  sienaFit$effects[-escreen,]$shortName[i])
+									  sienaFit$effects[-escreen,]$shortName[i] &
+									  dummyProps$interaction1 ==
+									  sienaFit$effects[-escreen,]$interaction1[i])
 						 if (length(tmp) > 0)
 						 {
 							doTests[tmp] <- TRUE
@@ -215,6 +227,7 @@
 							NA
 						 }
 					})
+
 	dim(effectTest) <- c(length(indBaseEffects), 1)
 	effectTestP <- round(1 - pchisq(effectTest, apply(toTest, 1, sum)), 5)
 	rownames(effectTestP) <- baseNames

Modified: pkg/RSienaTest/R/sienaeffects.r
===================================================================
--- pkg/RSienaTest/R/sienaeffects.r	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/R/sienaeffects.r	2010-05-29 21:10:40 UTC (rev 89)
@@ -11,7 +11,8 @@
 ##@includeEffect DataCreate
 includeEffects <- function(myeff, ..., include=TRUE, name=myeff$name[1],
                            type="eval", interaction1="", interaction2="",
-                           character=FALSE)
+                           character=FALSE, initialValue=0, fix=FALSE,
+						   timeDummy=",")
 {
 
     if (character)
@@ -40,6 +41,9 @@
     myeff$interaction1 == interaction1 &
     myeff$interaction2 == interaction2
     myeff[use, "include"] <- include
+	myeff[use, "fix"] <- fix
+	myeff[use, "initialValue"] <- initialValue
+	myeff[use, "timeDummy"] <- timeDummy
     print.data.frame(myeff[use, c("name", "shortName", "type", "interaction1",
                        "interaction2", "include")])
     myeff

Modified: pkg/RSienaTest/doc/RSiena.bib
===================================================================
--- pkg/RSienaTest/doc/RSiena.bib	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/doc/RSiena.bib	2010-05-29 21:10:40 UTC (rev 89)
@@ -704,7 +704,7 @@
 
 @ARTICLE{KiesnerEA04,
   AUTHOR="Jeff {Kiesner}, Margaret {Kerr}, H{\aa}kan {Stattin} ",
-  TITLE = "‘‘Very Important Persons'' in adolescence: going beyond
+  TITLE = "ï¿½ï¿½Very Important Persons'' in adolescence: going beyond
           in-school, single friendships in the study of peer homophily",
     YEAR = "2004",
   JOURNAL = "Journal of Adolescence",
@@ -894,6 +894,24 @@
   volume = 	 20,
   pages = 	 {1--21}}
 
+ at ARTICLE{Lospinoso2010a,
+	TITLE={Assessing and Accounting for Time Heterogeneity in Stochastic Actor Oriented Models},
+	AUTHOR={Joshua A. Lospinoso and M. Schweinberger and T.A.B. Snijders and R.M. Ripley},
+	JOURNAL={Advances in Data Analysis and Computation},
+	VOLUME={Special Issue on Social Networks (Submitted)},
+	YEAR={2010},
+	URL = {http://www.stats.ox.ac.uk/~lospinos}
+}
+ at Article{Lospinoso2010b,
+	TITLE = {Testing and Modeling Time Heterogeneity in Longitudinal Studies of Social Networks: A Tutorial in RSiena},
+	AUTHOR = {Joshua A. Lospinoso}, 
+	JOURNAL = {Connections},
+	VOLUME = {In progress.},
+	YEAR = {2010},
+	URL = {http://www.stats.ox.ac.uk/~lospinos}
+}
+
+
 @InCollection{Lepkowski89,
   author = 		 {J. M. Lepkowski},
   title = 		 {Treatment of wave nonresponse in panel surveys.},
@@ -913,6 +931,9 @@
 EDITOR = "Philipp {Herder-Dorneich} and Karl-Ernst {Schenk} and Dieter {Schmidtchen}",
 ADDRESS = "T{\"{u}}bingen"
 }
+
+
+
  @INCOLLECTION{MacyEA03,
   AUTHOR = "M.W. {Macy} and J. {Kitts} and A. {Flache} and S. {Benard}",
   TITLE = "A {H}opfield model of emergent structure",
@@ -1582,12 +1603,12 @@
   AUTHOR="Frans N. {Stokman}",
   TITLE = "Was verbindet uns wann mit wem? Inhalt und Struktur in der Analyse sozialer Netzwerke",
     YEAR = "2005",
-  JOURNAL = "Kölner Zeitschrift für Soziologie"
+  JOURNAL = "Kï¿½lner Zeitschrift fï¿½r Soziologie"
 }
 
  @INCOLLECTION{StokmanZeggelink96,
   AUTHOR = "Frans N. {Stokman}, Evelien P.H. {Zeggelink}",
-  TITLE = "‘Self-organizing' friendship networks",
+  TITLE = "ï¿½Self-organizing' friendship networks",
   BOOKTITLE = "Frontiers in social dilemmas research",
   PUBLISHER = "Springer",
   YEAR = "1996",
@@ -1620,7 +1641,7 @@
    Author = {ter Braak, C.J.F.},
    Title = {Permutation versus Bootstrap Significance Tests in Mulitple Regression and ANOVA.},
    BookTitle = {Bootstrapping and Related Techniques},
-   Editor = {Jöckel, K-H. and Rothe, G. and Sendler, W.},
+   Editor = {Jï¿½ckel, K-H. and Rothe, G. and Sendler, W.},
    Publisher = {Springer-Verlag},
    Address = {Berlin},
    Year = {1992} }

Modified: pkg/RSienaTest/doc/s_man400.tex
===================================================================
--- pkg/RSienaTest/doc/s_man400.tex	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/doc/s_man400.tex	2010-05-29 21:10:40 UTC (rev 89)
@@ -767,7 +767,8 @@
 #####################################GENERAL###################################
 
 # This is an R script for getting started with RSiena, written by
-# Robin Gauthier, Tom Snijders, and Ruth Ripley.
+# Robin Gauthier, Tom Snijders, and Ruth Ripley, with contributions
+# from Josh Lospinoso.
 # Lines starting with # are not processed by R but treated as comments.
 
 # R is case sensitive.
@@ -1122,6 +1123,44 @@
 # previous run from which to extract the theta estimates and the
 derivatives. Phase 1 will be omitted in this case.
 
+# Suppose also that you have a network with more than two time periods.
+# There are facilities available for dealing with such situations, where
+# effects may differ over time.
+#
+# Load a dataset with three time periods:
+
+	mynet2 <- sienaNet(array(c(s501, s502, s503), dim=c(50, 50, 3)))
+
+# Set up your data objects and effects objects as before:
+
+	mydata2 <- sienaDataCreate(mynet2)
+	myeff2 <- getEffects(mydata2)
+
+# Try out the new includeEffects functionality for adding effects:
+
+	myeff2 <- includeEffects(myeff2, transTrip, balance)
+
+# Run an estimation with your old model parameters:
+
+	mymodel2 <- sienaModelCreate(fn=simstats0c, nsub=4, n3=1000)
+	ans <- siena07(mymodel2, data=mydata2, effects=myeff2, batch=TRUE)
+
+# Conduct the score type test to assess whether heterogeneity is present.
+
+	tt <- sienaTimeTest(ans)
+
+# Try a plot to see a visual representation of the time test:
+
+	plot(tt, effects=1:4, dims=c(2,2))
+
+# Add a few time dummy terms:
+
+	myeff3<- includeTimeDummy(myeff2, recip, balance, timeDummy="2")
+
+# And estimate the new model:
+
+ans3 <- siena07(mymodel2, data=mydata2, effects=myeff2, batch=TRUE)
+
 #######################VIEWING THE NETWORK IN R#######################
 
 # We can make connections with other R packages, e.g., Carter Butts's
@@ -1349,7 +1388,9 @@
 
 More specific references are \citet{Schweinberger10} for the score-type goodness
 of fit tests and \citet{SchweinbergerSnijders07a} for the calculation of
-standard errors of the Method of Moments estimators .
+standard errors of the Method of Moments estimators. For assessing and
+correcting time heterogeneity, and associated model selection considerations,
+refer to \citet*{Lospinoso2010a, Lospinoso2010b}.
 
 A tutorial is \citet*{SnijdersEA10b}.
 
@@ -3321,6 +3362,32 @@
 %of the unobserved actor heterogeneity.
 \fi
 
+\subsection{Time heterogeneity in model parameters}
+Currently for the case of a one mode network, you can include
+time heterogeneous parameters in your model. Consider the reformulation of
+the evaluation function into
+\begin{align}
+f^{(m)}_{ij}(\mathbf{x})= \sum_k (\beta_k + \delta_k^{(m)} h_k^{(m)}) s_{ik}(\mathbf{x}(i \leadsto j))
+\label{eq:fmij}
+\end{align}
+where $\delta_k^{(m)}$ are the time dummy interacted effect parameters. You
+can include these in your model simply via the function
+\begin{verbatim}
+myeffects <- includeTimeDummy(myeffects, density, reciprocity, timeDummy="2,3,6")
+\end{verbatim}
+which would add three time dummy terms to each effect listed in the function.
+
+We recommend that you start with simple models, and use the score type test for
+assessing heterogeneity, i.e.
+\begin{verbatim}
+tt <- sienaTimeTest(ans)
+\end{verbatim}
+to decide which dummy terms to include.
+
+See \citet{Lospinoso2010a} for a technical presentation of how the test works,
+and \citet{Lospinoso2010b} for a walkthrough on model selection.
+
+
 \begin{print}
 \newpage
 \end{print}
@@ -4290,8 +4357,61 @@
 which under the null hypothesis of equal parameters has an approximating
 standard normal distribution.
 
+\newpage
+\subsection{Testing time heterogeneity in parameters}
+We initially assume that $\beta$ does not vary over time, yielding a \emph{restricted model}. Our data contains $|\mathcal{M}|$ observations, and we estimate the restricted model the method of moments. We wish to test whether the \emph{restricted model} is misspecified with respect to time heterogeneity. Formally, define a vector of time dummy terms $\mathbf{h}$:
+\begin{align}
+h_k^{(m)}=\left\{
+\begin{array}{ll}
+1& \{m : w_m \in \mathcal{W}, m \neq 1\}\\
+0& \mbox{~elsewhere~}
+\end{array}
+\right . ,
+\end{align}
+where $k$ corresponds to an effect included in the model.\footnote{The dummy $\delta_k^{(1)}$ is always zero so that period $w_1$ is (arbitrarily) considered the reference period.} An \emph{unrestricted model} which allows for time heterogeneity in all of the effects is considered as a modification of \eqref{eq:fij}:
 
+\begin{align}
+f^{(m)}_{ij}(\mathbf{x})= \sum_k (\beta_k + \delta_k^{(m)} h_k^{(m)}) s_{ik}(\mathbf{x}(i \leadsto j))
+\label{eq:fmij}
+\end{align}
+where $\delta_k^{(m)}$ are the time dummy interacted effect parameters. One way to formulate the testing problem of assessing time hetergeneity is the following:
+\begin{align}
+H_0:\delta_k^{(m)} = 0 & \mbox{~for all~} k,m \notag\\
+H_1:\delta_k^{(m)} \neq 0 & \mbox{ for some } k,m .
+\label{hyptest}
+\end{align}
 
+An application of the score test is given for the special case of parameter heterogeneity
+by \citet{Lospinoso2010a} and implemented in RSiena.
+To apply the test to your dataset, run an estimation in the usual way, e.g.:
+\begin{verbatim}
+mymodel <- sienaModelCreate(fn=simstats0c, nsub=2, n3=100)
+mynet1 <- sienaNet(array(c(s501, s502, s503), dim=c(50, 50, 3)))
+mydata <- sienaDataCreate(mynet1)
+myeff <- getEffects(mydata)
+myeff <- includeEffects(myeff, transTrip, balance)
+ans <- siena07(mymodel, data=mydata, effects=myeff, batch=TRUE)
+\end{verbatim}
+and conduct the timetest through
+\begin{verbatim}
+## Conduct the score type test to assess whether heterogeneity is present.
+tt <- sienaTimeTest(ans)
+plot(tt, effects=1:2, dim=c(1,2))
+\end{verbatim}
+If it is decided that you wish to add time dummy terms, you may do so via
+\begin{verbatim}
+myeff <- includeTimeDummy(myeff, recip, balance, timeDummy="2")
+ans2 <- siena07(mymodel, data=mydata, effects=myeff, batch=TRUE)
+\end{verbatim}
+and test again,
+\begin{verbatim}
+tt2 <- sienaTimeTest(ans2)
+\end{verbatim}
+and so on.
+
+See \citet{Lospinoso2010b} for a walkthrough of the model selection process
+for time dummy terms.
+
 \begin{print}
 \newpage
 \end{print}
@@ -7259,8 +7379,13 @@
 (Programmers should consult the changeLog file on CRAN or in the R-forge
 repository.)
 \begin{itemize}
+\item 2010-04-24 R-forge revision 82-88
+New functionality for sienaTimeTest, new function includeTimeDummy, added funcitonality
+to includeEffects, sienaDataCreate.
 \item 2010-04-24 R-forge revision 81
 New print, summary and edit methods for Siena effects objects
+\item 2010-04-24 R-forge revision 81
+New print, summary and edit methods for Siena effects objects
 \item 2010-04-24 R-forge revision 80
 \begin{itemize}
 \item fixed bug causing crash with rate effects and bipartite networks.

Modified: pkg/RSienaTest/man/includeEffects.Rd
===================================================================
--- pkg/RSienaTest/man/includeEffects.Rd	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/man/includeEffects.Rd	2010-05-29 21:10:40 UTC (rev 89)
@@ -9,7 +9,8 @@
 }
 \usage{
 includeEffects(myeff, ..., include = TRUE, name = myeff$name[1],
- type = "eval", interaction1 = "", interaction2 = "", character=FALSE)
+ type = "eval", interaction1 = "", interaction2 = "", character=FALSE,
+ initialValue=0, fix=FALSE, timeDummy=",")
 }
 \arguments{
   \item{myeff}{a Siena effects object as created by \code{\link{getEffects}}
@@ -35,6 +36,10 @@
   effects e.g. covariate name or behavior variable name.
 }
 \item{character}{Boolean: are the effect names character strings or not}
+\item{initialValue}{Float: Set an initial value for the effect.}
+\item{fix}{Boolean: Fix the effect during estimation}
+\item{timeDummy}{string: Comma delimited string of which periods to
+dummy. Alternately, use \code{\link{includeTimeDummy}}.}
 }
 \details{The arguments should identify the effects completely. The
   include column will be set to the value requested (TRUE or FALSE).

Modified: pkg/RSienaTest/man/sienaTimeTest.Rd
===================================================================
--- pkg/RSienaTest/man/sienaTimeTest.Rd	2010-05-28 17:27:56 UTC (rev 88)
+++ pkg/RSienaTest/man/sienaTimeTest.Rd	2010-05-29 21:10:40 UTC (rev 89)
@@ -1,5 +1,6 @@
 \name{sienaTimeTest}
 \alias{sienaTimeTest}
+\alias{includeTimeDummy}
 \alias{plot.sienaTimeTest}
 \title{Functions to assess and account for time heterogeneity of parameters}
 \description{
@@ -25,7 +26,7 @@
  }
 
 \usage{
-timetest <- sienaTimeTest(sienaFit)
+timetest <- sienaTimeTest(sienaFit, condition=FALSE)
 \method{plot}{sienaTimeTest}(x, pairwise=FALSE, effects=1:2,
 	dims=c(2, 1), scale=.2, plevels=c(.1, .05, .025), multiplot=FALSE, ...)
 myeffects <- includeTimeDummy(sienaEffects, outdegree,
@@ -33,6 +34,8 @@
 }
 \arguments{
   \item{sienaFit}{ A sienaFit object returned by \code{siena07}. }
+  \item{condition}{ Whether to orthogonalize individual score tests
+  against base effects and un-estimated dummy terms or just base effects.}
   \item{x}{ A sienaTimeTest object returned by \code{sienaTimeTest}. }
   \item{pairwise}{ A Boolean value corresponding to whether the user
   would like a pairwise plot of the simulated statistics to assess
@@ -56,7 +59,6 @@
  \item{\dots}{For extra arguments (none used at present)}
 }
 \details{
-
   This test follows the score type test of Schweinburger (2007b) as
  implemented by Lospinoso et. al. (2010) by using statistics already
  calculated at each wave to populate vectors of partitioned moment
@@ -98,7 +100,6 @@
  \code{timeDummy} column provides facilities within a sienaFit object to
  quickly adjust the model and re-estimate in the usual way. Alternately,
  the includeTimeDummy function can be used to add the desired dummies.
-
 }
 \value{
   \code{sienaTimeTest} Returns a list containing many items,
@@ -110,6 +111,7 @@
 	step estimates, and a p-value for \code{H0:} the unrestricted
 	parameters are equal to zero.}
 }
+
 \references{See \url{http://www.stats.ox.ac.uk/~snijders/siena/}
   for general information on RSiena.
 
@@ -124,9 +126,10 @@
 51(9):4465-4483, 2007.  ISSN 0167-9473.
 
 M. Schweinberger. Chapter 4: Statistical Modeling of Network Panel Data:
-Goodness of Fit. PhD thesis, University of Groningen, 2007.}
+Goodness of Fit. PhD thesis, University of Groningen, 2007.
+}
 \author{Josh Lospinoso}
-\seealso{\code{\link{siena07}}, \code{\link{sienaTimeTest}}}
+\seealso{\code{\link{siena07}}, \code{\link{sienaTimeTest}}, \code{\link{includeTimeDummy}}}
 \examples{
 ## Estimate a restricted model
 mymodel <- sienaModelCreate(fn=simstats0c, nsub=2, n3=100)
@@ -184,4 +187,5 @@
 			interaction1="myccov")
 ans <- siena07(mymodel, data=mydata, effects=myeff, batch=TRUE)
 }
+}
 \keyword{models}

