[Rsiena-commits] r5 - in pkg/RSiena: . R inst/doc man src src/data src/model src/model/effects src/model/tables
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Jul 30 23:15:50 CEST 2009
Author: ripleyrm
Date: 2009-07-30 23:15:49 +0200 (Thu, 30 Jul 2009)
New Revision: 5
Modified:
pkg/RSiena/DESCRIPTION
pkg/RSiena/R/effectsInfo.R
pkg/RSiena/R/interrupt.r
pkg/RSiena/R/phase1.r
pkg/RSiena/R/phase2.r
pkg/RSiena/R/phase3.r
pkg/RSiena/R/robmon.r
pkg/RSiena/R/siena07.r
pkg/RSiena/R/sienaDataCreate.r
pkg/RSiena/R/simstatsc.r
pkg/RSiena/inst/doc/s_man400.pdf
pkg/RSiena/man/RSiena-package.Rd
pkg/RSiena/man/model.create.Rd
pkg/RSiena/man/siena07.Rd
pkg/RSiena/man/simstats0c.Rd
pkg/RSiena/src/data/ChangingDyadicCovariate.cpp
pkg/RSiena/src/data/ChangingDyadicCovariate.h
pkg/RSiena/src/data/ConstantDyadicCovariate.cpp
pkg/RSiena/src/data/ConstantDyadicCovariate.h
pkg/RSiena/src/data/IncidentTieIterator.cpp
pkg/RSiena/src/data/IncidentTieIterator.h
pkg/RSiena/src/data/Network.cpp
pkg/RSiena/src/data/Network.h
pkg/RSiena/src/data/OneModeNetwork.cpp
pkg/RSiena/src/data/OneModeNetwork.h
pkg/RSiena/src/data/OneModeNetworkLongitudinalData.cpp
pkg/RSiena/src/model/EpochSimulation.cpp
pkg/RSiena/src/model/effects/AllEffects.h
pkg/RSiena/src/model/effects/CovariateAlterEffect.cpp
pkg/RSiena/src/model/effects/CovariateAlterEffect.h
pkg/RSiena/src/model/effects/CovariateEgoAlterEffect.cpp
pkg/RSiena/src/model/effects/CovariateEgoAlterEffect.h
pkg/RSiena/src/model/effects/CovariateSimilarityEffect.cpp
pkg/RSiena/src/model/effects/CovariateSimilarityEffect.h
pkg/RSiena/src/model/effects/DyadicCovariateDependentNetworkEffect.cpp
pkg/RSiena/src/model/effects/DyadicCovariateDependentNetworkEffect.h
pkg/RSiena/src/model/effects/EffectFactory.cpp
pkg/RSiena/src/model/effects/NetworkEffect.cpp
pkg/RSiena/src/model/effects/SameCovariateEffect.cpp
pkg/RSiena/src/model/effects/SameCovariateEffect.h
pkg/RSiena/src/model/tables/TwoPathTable.cpp
pkg/RSiena/src/siena07.cpp
Log:
New effects, bug fixes
Modified: pkg/RSiena/DESCRIPTION
===================================================================
--- pkg/RSiena/DESCRIPTION 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/DESCRIPTION 2009-07-30 21:15:49 UTC (rev 5)
@@ -1,8 +1,8 @@
Package: RSiena
Type: Package
Title: Siena - Simulation Investigation for Empirical Network Analysis
-Version: 1.0.2
-Date: 2009-07-8
+Version: 1.0.3
+Date: 2009-07-30
Author: Various
Depends: R (>= 2.7.0)
Imports: Matrix
@@ -14,4 +14,4 @@
LazyLoad: yes
LazyData: yes
URL: http://www.stats.ox.ac.uk/~snijders/siena
-Packaged: 2009-07-08 19:35:50 UTC; ruth
+Packaged: 2009-07-30 20:55:58 UTC; ruth
Modified: pkg/RSiena/R/effectsInfo.R
===================================================================
--- pkg/RSiena/R/effectsInfo.R 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/effectsInfo.R 2009-07-30 21:15:49 UTC (rev 5)
@@ -47,8 +47,8 @@
"nbrDist2", "nbrDist2twice", "denseTriads", "inPop", "inPopSqrt",
"outPop", "outPopSqrt", "inAct", "inActSqrt", "outAct", "outActSqrt",
"outInv", "outSqInv", "outOutAss", "outInAss", "inOutAss", "inInAss"
-), parm = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
-0, 0, 0, 0, 2, 2, 2, 2)), .Names = c("EffectName", "FunctionName",
+), parm = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 5, 0, 0, 0, 0, 0, 0,
+0, 0, 1, 1, 2, 2, 2, 2)), .Names = c("EffectName", "FunctionName",
"Endowment.", "ShortName", "parm"), row.names = c(NA, -25L), class = "data.frame")
symmetricObjEffects <-
structure(list(EffectName = c("degree (density)", "transitive triads",
Modified: pkg/RSiena/R/interrupt.r
===================================================================
--- pkg/RSiena/R/interrupt.r 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/interrupt.r 2009-07-30 21:15:49 UTC (rev 5)
@@ -25,69 +25,69 @@
if (is.null(tt))
{
library(tcltk)
- tt<- tktoplevel()
+ tt <- tktoplevel()
}
tkwm.title(tt,'Siena07')
- frame<- tkframe(tt,width=300,height=300,relief='ridge')
- tkpack(frame,side='top',padx=5)
- button1<- tkbutton(frame,command=myInterrupt,text='Interrupt')
- button2<- tkbutton(frame,command=myEndPhase2,text='End Phase2',
+ frame <- tkframe(tt, width=300, height=300, relief='ridge')
+ tkpack(frame, side='top', padx=5)
+ button1 <- tkbutton(frame, command=myInterrupt, text='Interrupt')
+ button2 <- tkbutton(frame, command=myEndPhase2, text='End Phase2',
state='disabled')
- button3<- tkbutton(frame,command=myRestart,text='Restart')
- tkgrid.configure(button1,column=1,columnspan=2,row=1,padx=20,pady=20)
- tkgrid.configure(button2,column=3,row=1,padx=20)
- tkgrid.configure(button3,column=4,row=1,padx=20)
- phaselabel<-tklabel(frame,text='Phase')
- subphaselabel<-tklabel(frame,text='Subphase',state='disabled')
- iterationlabel<-tklabel(frame,text='Iteration')
- label1<-tklabel(frame,text='ProgressBar')
+ button3 <- tkbutton(frame, command=myRestart, text='Restart')
+ tkgrid.configure(button1, column=1, columnspan=2, row=1, padx=20, pady=20)
+ tkgrid.configure(button2, column=3, row=1, padx=20)
+ tkgrid.configure(button3, column=4, row=1, padx=20)
+ phaselabel <- tklabel(frame, text='Phase')
+ subphaselabel <- tklabel(frame, text='Subphase', state='disabled')
+ iterationlabel <- tklabel(frame, text='Iteration')
+ label1 <- tklabel(frame, text='ProgressBar')
- phase<- tkentry(frame, width=2)
+ phase <- tkentry(frame, width=2)
- subphase<- tkentry(frame, width=2,state='disabled')
- iteration<- tkentry(frame, width=6)
- progressbar<- ttkprogressbar(frame,max=2000,length=120)
+ subphase <- tkentry(frame, width=2, state='disabled')
+ iteration <- tkentry(frame, width=6)
+ progressbar <- ttkprogressbar(frame, max=2000, length=120)
- tkgrid.configure(phaselabel,column=1,row=2,pady=5)
- tkgrid.configure(subphaselabel,column=2,row=2)
- tkgrid.configure(iterationlabel,column=3,row=2)
- tkgrid.configure(label1,column=4,row=2,padx=5)
- tkgrid.configure(phase,column=1,row=3,pady=3)
- tkgrid.configure(subphase,column=2,row=3,padx=10)
- tkgrid.configure(iteration,column=3,row=3,padx=10)
+ tkgrid.configure(phaselabel, column=1, row=2, pady=5)
+ tkgrid.configure(subphaselabel, column=2, row=2)
+ tkgrid.configure(iterationlabel, column=3, row=2)
+ tkgrid.configure(label1, column=4, row=2, padx=5)
+ tkgrid.configure(phase, column=1, row=3, pady=3)
+ tkgrid.configure(subphase, column=2, row=3, padx=10)
+ tkgrid.configure(iteration, column=3, row=3, padx=10)
- tkgrid.configure(progressbar,column=4,padx=5,row=3)
- label2<- tklabel(frame,text='Current parameter values')
- label3<- tklabel(frame,text='Quasi-autocorrelations')
- label4<- tklabel(frame,text='Deviation values')
+ tkgrid.configure(progressbar, column=4, padx=5, row=3)
+ label2 <- tklabel(frame, text='Current parameter values')
+ label3 <- tklabel(frame, text='Quasi-autocorrelations')
+ label4 <- tklabel(frame, text='Deviation values')
- tkgrid.configure(label2,column=1,columnspan=2,row=4,padx=10)
- tkgrid.configure(label3,column=3,row=4,padx=10)
- tkgrid.configure(label4,column=4,row=4,padx=10)
+ tkgrid.configure(label2, column=1, columnspan=2, row=4, padx=10)
+ tkgrid.configure(label3, column=3, row=4, padx=10)
+ tkgrid.configure(label4, column=4, row=4, padx=10)
- text1<- tktext(frame, height=6,width=14)
+ text1 <- tktext(frame, height=6, width=14)
- text2<- tktext(frame, height=6,width=14)
- text3<- tktext(frame, height=6,width=14)
- tkgrid.configure(text1,column=1,columnspan=2,row=5,padx=20,pady=5)
+ text2 <- tktext(frame, height=6, width=14)
+ text3 <- tktext(frame, height=6, width=14)
+ tkgrid.configure(text1, column=1, columnspan=2, row=5, padx=20, pady=5)
- tkgrid.configure(text2,column=3,row=5,padx=20)
- tkgrid.configure(text3,column=4,row=5,padx=20)
+ tkgrid.configure(text2, column=3, row=5, padx=20)
+ tkgrid.configure(text3, column=4, row=5, padx=20)
ilcampo <- tclVar()
- tcl("image","create","photo",ilcampo,file=imagepath)
- frame2<- tkframe(tt,width=300,height=300,relief='ridge')
- tkpack(frame2,side='bottom',padx=5)
- imgAsLabel <- tklabel(frame2,image=ilcampo)
- tkgrid.configure(imgAsLabel,pady=10)
- tkinsert(phase,0,' 1')
+ tcl("image", "create", "photo", ilcampo, file=imagepath)
+ frame2 <- tkframe(tt, width=300, height=300, relief='ridge')
+ tkpack(frame2, side='bottom', padx=5)
+ imgAsLabel <- tklabel(frame2, image=ilcampo)
+ tkgrid.configure(imgAsLabel, pady=10)
+ tkinsert(phase, 0, ' 1')
tkgrab.set(tt)
tcl('update')
# browser()
tkfocus(tt)
# cat('here\n')
- list(tt=tt,pb=progressbar,earlyEndPhase2=button2,current=text1,
+ list(tt=tt, pb=progressbar, earlyEndPhase2=button2, current=text1,
quasi=text2, deviations=text3, phase=phase, subphase=subphase,
- iteration=iteration,subphaselabel=subphaselabel)
+ iteration=iteration, subphaselabel=subphaselabel)
}
Modified: pkg/RSiena/R/phase1.r
===================================================================
--- pkg/RSiena/R/phase1.r 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/phase1.r 2009-07-30 21:15:49 UTC (rev 5)
@@ -14,14 +14,34 @@
## *****************************************************************************/
##args: x model object (readonly), z control object
##
-phase1.1<- function(z,x,...)
+phase1.1 <- function(z, x, ...)
{
## initialise phase 1
- z$SomeFixed<- FALSE
+ z$SomeFixed <- FALSE
z$Phase <- 1
- z<-AnnouncePhase(z,x)
+ f <- FRANstore()
+ int <- z$int
+ z <- AnnouncePhase(z, x)
+ z$phase1Its <- 0
+ ## fix up iteration numbers if using multiple processors
+ if (10 %% int == 0)
+ {
+ firstNit <- 10
+ }
+ else
+ {
+ firstNit <- 10 + int - 10 %% int
+ }
+ if ((z$n1 - firstNit) %% int == 0)
+ {
+ endNit <- z$n1
+ }
+ else
+ {
+ endNit <- z$n1 + int - (z$n1 - firstNit) %% int
+ }
+ z$n1 <- endNit
z$sf <- matrix(0, nrow = z$n1, ncol = z$pp)
- f <- FRANstore() ## the derivatives are model dependent!
z$sf2 <- array(0, dim=c(z$n1, f$observations - 1, z$pp))
z$ssc <- array(0, dim=c(z$n1, f$observations - 1, z$pp))
z$sdf <- array(0, dim=c(z$n1, z$pp, z$pp))
@@ -29,44 +49,56 @@
z$writefreq <- 10
z$DerivativeProblem <- FALSE
z$Deriv <- !z$FinDiff.method ## can do both in phase 3 but not here!
- for (nit in 1 : 10)
+ zsmall <- NULL
+ zsmall$theta <- z$theta
+ zsmall$Deriv <- z$Deriv
+ zsmall$Phase <- z$Phase
+ zsmall$nit <- z$nit
+ xsmall<- NULL
+ xsmall$cconditional <- x$cconditional
+ zsmall$condvar <- z$condvar
+ nits <- seq(1, firstNit, int)
+ nits6 <- min(nits[nits >= 6 ])
+ for (nit in nits)
{
- z$nit <- nit
- if (nit == 2)
- {
- time1 <- proc.time()['elapsed']
- if (x$checktime)
- {
- z$ctime <- time1
- }
- }
- if (nit==6)
- {
- time1<- proc.time()['elapsed']-time1
- if (time1>1e-5)
- {
- z$writefreq <- round(20.0/time1)
- }
- else
- {
- z$writefreq <- 5
- }
- if (is.batch())
- {
- z$writefreq <- 10*z$writefreq
- }
- z$writefreq<-roundfreq(z$writefreq)
- }
- #########################################
- ### do iteration
- #########################################
- z<-doPhase1it(z,x,...)
- #########################################
- ###
- #########################################
- if (!z$OK || UserRestartFlag() || UserInterruptFlag())
- return(z)
- }
+ z$nit <- nit
+ if (nit == nits[2])
+ {
+ time1 <- proc.time()['elapsed']
+ if (x$checktime)
+ {
+ z$ctime <- time1
+ }
+ }
+ if (nit == nits6)
+ {
+ time1 <- proc.time()['elapsed']-time1
+ if (time1 > 1e-5)
+ {
+ z$writefreq <- round(20.0/time1)
+ }
+ else
+ {
+ z$writefreq <- 5
+ }
+ if (is.batch())
+ {
+ z$writefreq <- 10 * z$writefreq
+ }
+ z$writefreq <- roundfreq(z$writefreq)
+ }
+ ## #######################################
+ ## # do iteration
+ ## #######################################
+ z <- doPhase1it(z, x, cl=z$cl, int=z$int, zsmall=zsmall,
+ xsmall=xsmall, ...)
+ ## #######################################
+ ## #
+ ## #######################################
+ if (!z$OK || UserRestartFlag() || UserInterruptFlag())
+ return(z)
+ }
+
if (z$FinDiff.method)
{
npos <- z$npos
@@ -75,42 +107,44 @@
j<- (1 : z$pp)[!z$fixed & npos < 5]
for (i in 1 : length(j))
{
- Report(c('After 10 iterations, only', npos[j[i]],
- 'differences in coordinate', j[i], '.\n'), cf)
- Report(c('with epsilon = ',z$epsilon[j[i]], '.\n'), cf)
+ Report(c("After ", z$n, " iterations, only ", npos[j[i]],
+ " differences in coordinate ", j[i], ".\n"), sep="",
+ cf)
+ Report(c("with epsilon = ", z$epsilon[j[i]], ".\n"), sep="", cf)
}
- use<- !z$fixed & npos <= 2
- z$epsilon[use]<- ifelse(z$posj[use],
- 3.0 * z$epsilon[use],
- 10.0 * z$epsilon[use])
- use<- !z$fixed & npos > 2 & npos < 5
- z$epsilon[use]<- ifelse(z$posj[use],
+ use <- !z$fixed & npos <= 2
+ z$epsilon[use] <- ifelse(z$posj[use],
+ 3.0 * z$epsilon[use],
+ 10.0 * z$epsilon[use])
+ use <- !z$fixed & npos > 2 & npos < 5
+ z$epsilon[use] <- ifelse(z$posj[use],
2.0 * z$epsilon[use],
sqrt(10.0) * z$epsilon[use])
- use<- !z$fixed & npos < 5
+ use <- !z$fixed & npos < 5
z$epsilon[use] <- pmin(100.0 * z$scale[use], z$epsilon[use])
z$epsilon[use] <- pmax(0.1 * z$scale[use], z$epsilon[use])
- Report(c('New epsilon =', z$epsilon[use],'.\n'),cf,fill=80)
+ Report(c("New epsilon =", paste(" ", z$epsilon[use], collapse=""),
+ "\n"), sep="", cf, fill=80)
if (z$repeatsForEpsilon <= 4)
{
- Report('Change value of epsilon and restart Phase 1.\n',cf)
- z$epsilonProblem<- TRUE
+ Report("Change value of epsilon and restart Phase 1.\n", cf)
+ z$epsilonProblem <- TRUE
return(z)
}
- if (any(npos[!z$fixed]<=1))
+ if (any(npos[!z$fixed] <= 1))
{
j0s <- which(!z$fixed & npos <= 1)
for (j0 in j0s)
{
- Report (c('Difficulties with parameter', j0, '.\n'),
+ Report (c("Difficulties with parameter", j0, ".\n"),
outf)
- Report(c('After 10 iterations, only', npos[j0],
- 'differences in this coordinate.\n'), outf)
- Report('This parameter is fixed and not estimated.\n',
+ Report(c("After ", z$n, " iterations, only", npos[j0],
+ "differences in this coordinate.\n"), outf)
+ Report("This parameter is fixed and not estimated.\n",
outf)
- Report(c('Fix parameter',j0,',in phase 1.\n'), cf)
+ Report(c("Fix parameter", j0, ", in phase 1.\n"), cf)
}
- z$newFixed[j0s]<- TRUE
+ z$newFixed[j0s] <- TRUE
z$fixed[j0s] <- TRUE
z$SomeFixed <- TRUE
}
@@ -120,39 +154,73 @@
}
}
}
-
z
}
-doPhase1it<- function(z, x, ...)
+doPhase1it<- function(z, x, cl, int, zsmall, xsmall, ...)
{
- z$n <- z$n + 1
DisplayIteration(z)
- ## seed0<- .Random.seed ##store seed for use in finitedifferences
- ## now stored on z as an array by period and reused in finite differences
- zz <- x$FRAN(z, x, INIT=FALSE, ...)
- fra <- colSums(zz$fra)
- if (!zz$OK)
+ if (int == 1)
{
- z$OK <- zz$OK
- return(z)
+ zz <- x$FRAN(zsmall, xsmall)
+ if (!zz$OK)
+ {
+ z$OK <- zz$OK
+ z$zz <- zz
+ return(z)
+ }
+ z$n <- z$n + 1
+ z$phase1Its <- z$phase1Its + 1
}
- if (!is.null(zz[['sc']]))
+ else
{
- z$ssc[z$nit, , ] <- zz$sc
+ zz <- clusterCall(cl, usesim, zsmall, xsmall)
+ z$n <- z$n + z$int
+ z$phase1Its <- z$phase1Its + int
}
- fra <- fra -z$targets
- z$sf[z$nit, ] <- fra
- z$sf2[z$nit, , ] <- zz$fra
+ if (int == 1)
+ {
+ fra <- colSums(zz$fra)
+ fra <- fra - z$targets
+ z$sf[z$nit, ] <- fra
+ z$sf2[z$nit, , ] <- zz$fra
+ }
+ else
+ {
+ for (i in 1:int)
+ {
+ fra <- colSums(zz[[i]]$fra)
+ fra <- fra - z$targets
+ z$sf[z$nit + (i - 1), ] <- fra
+ z$sf2[z$nit + (i - 1), , ] <- zz[[i]]$fra
+ }
+
+ }
if (z$FinDiff.method)
{
- z <- FiniteDifferences(z, x, fra + z$targets, ...)
+ z <- FiniteDifferences(z, x, fra + z$targets, z$cl, z$int, ...)
z$sdf[z$nit, , ] <- z$sdf0
}
else if (x$maxlike)
{
z$sdf[z$nit, , ] <- zz$dff
}
+ else
+ {
+ if (int==1)
+ {
+ if (!is.null(zz[['sc']]))
+ z$ssc[z$nit , ,] <- zz$sc
+ }
+ else
+ {
+ for (i in 1:int)
+ {
+ if (!is.null(zz[[i]][['sc']]))
+ z$ssc[z$nit + (i - 1), , ] <- zz[[i]]$sc
+ }
+ }
+ }
CheckBreaks()
if (UserInterruptFlag() || UserRestartFlag())
{
@@ -171,7 +239,8 @@
z$pb <- setProgressBar(z$pb, val)
progress <- val / z$pb$pbmax * 100
if (z$nit <= 5 || z$nit %% z$writefreq == 0 || z$nit %%5 == 0 ||
- x$maxlike || x$FinDiff.method)
+ x$maxlike || x$FinDiff.method ||
+ (int > 1 && z$nit %% z$writefreq < int))
{
# Report(c('Phase', z$Phase, 'Iteration ', z$nit, '\n'))
if (is.batch())
@@ -186,15 +255,26 @@
#browser()
z
}
-phase1.2<- function(z, x, ...)
+phase1.2 <- function(z, x, ...)
{
##finish phase 1 iterations and do end-of-phase processing
- if (z$n1 > 10)
+ zsmall <- NULL
+ zsmall$theta <- z$theta
+ zsmall$Deriv <- z$Deriv
+ zsmall$Phase <- z$Phase
+ zsmall$nit <- z$nit
+ xsmall<- NULL
+ xsmall$cconditional <- x$cconditional
+ zsmall$condvar <- z$condvar
+ int <- z$int
+ if (z$n1 > z$phase1Its)
{
- for (nit in 11 : z$n1)
+ nits <- seq((z$phase1Its+1), z$n1, int)
+ for (nit in nits)
{
z$nit <- nit
- z <- doPhase1it(z, x, ...)
+ z <- doPhase1it(z, x, cl=z$cl, int=int, zsmall=zsmall,
+ xsmall=xsmall, ...)
if (!z$OK || UserInterruptFlag() || UserRestartFlag())
{
return(z)
@@ -294,9 +374,10 @@
{
z$theta[!z$fixed] <- z$theta[!z$fixed] - fchange[!z$fixed]
}
- Report(c('Phase 1 achieved after ', z$nit, ' iterations.\n'), cf)
+ Report(c('Phase 1 achieved after ', z$phase1Its, ' iterations.\n'), cf)
WriteOutTheta(z)
- z$nitPhase1 <- z$nit
+ z$nitPhase1 <- z$phase1Its
+ z$phase1devs <- z$sf
##browser()
z
}
@@ -388,7 +469,9 @@
FiniteDifferences <- function(z, x, fra, cl, int=1, ...)
{
fras <- array(0, dim = c(int, z$pp, z$pp))
- # browser()
+ xsmall<- NULL
+ xsmall$cconditional <- x$cconditional
+ # browser()
for (i in 1 : z$pp)
{
zdummy <- z[c('theta', 'Deriv')]
@@ -397,10 +480,10 @@
zdummy$theta[i] <- z$theta[i] + z$epsilon[i]
}
## assign('.Random.seed',randomseed,pos=1)
- if (!z$Phase == 3 || int == 1)
+ if (int == 1)
{
- zz <- x$FRAN(zdummy, list(x$cconditional), INIT=FALSE,
- fromFiniteDiff=TRUE, ...)
+ zz <- x$FRAN(zdummy, xsmall, INIT=FALSE,
+ fromFiniteDiff=TRUE)
if (!zz$OK)
{
z$OK <- zz$OK
@@ -409,10 +492,10 @@
}
else
{
- zz <- clusterCall(cl, x$FRAN, zdummy, list(x$cconditional),
- INIT=FALSE, fromFiniteDiff=TRUE, ...)
+ zz <- clusterCall(cl, x$FRAN, zdummy, xsmall,
+ INIT=FALSE, fromFiniteDiff=TRUE)
}
- if (!z$Phase == 3 || int == 1)
+ if (int == 1)
{
fras[1, i, ] <- colSums(zz$fra) - fra
}
Modified: pkg/RSiena/R/phase2.r
===================================================================
--- pkg/RSiena/R/phase2.r 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/phase2.r 2009-07-30 21:15:49 UTC (rev 5)
@@ -11,11 +11,44 @@
## ****************************************************************************/
## args: z: internal control object
## x: model object (readonly as not returned)
-phase2.1<- function(z,x,...)
+usesim <- function(...)
{
+ simstats0c(...)
+}
+storeinFRANstore <- function(...)
+{
+ FRANstore(...)
+}
+phase2.1<- function(z, x, useCluster, noClusters, initC,...)
+{
#initialise phase2
+ int <- 1
+ # cl <- NULL
+ f <- FRANstore()
+ # if (useCluster)
+ # {
+ # if (!is.null(x$simstats0c) && !x$simstats0c)
+ # {
+ # stop("Multiple processors only for simstats0c at present")
+ # }
+ # cl <- makeCluster(rep("localhost", noClusters), type = "SOCK",
+ # outfile = 'cluster.out')
+ # clusterSetupRNG(cl, seed = rep(1, 6))
+ # clusterCall(cl, library, 'RSiena', character.only = TRUE)
+ # clusterCall(cl, storeinFRANstore, f)
+ # ans <- clusterCall(cl, storeinFRANstore)
+ # int <- noClusters
+ # if (initC)
+ # {
+ # ## ans <- clusterCall(cl, x$FRAN, z, x, INIT=FALSE, initC = initC)
+ # ans <- clusterCall(cl, usesim, z, x,
+ # INIT=FALSE, initC = initC)
+ # }
+ # }
+ # z$cl <- cl
+ # z$int <- int
z$Phase <- 2
- z$writefreq<- 1
+ z$writefreq <- 1
if (!is.batch())
{
tkconfigure(z$tkvars$earlyEndPhase2,state='normal')
@@ -72,6 +105,7 @@
## do the iterations for this repeat of this subphase
## ##############################################
z <- doIterations(z, x, subphase, ...)
+ ## if (z$nit == 50) browser()
if (!z$OK || UserInterruptFlag() || UserRestartFlag() ||
EarlyEndPhase2Flag())
{
@@ -167,10 +201,19 @@
{
z$nit <- 0
ac <- 0
+ zsmall <- NULL
+ zsmall$theta <- z$theta
+ zsmall$Deriv <- z$Deriv
+ zsmall$Phase<- z$Phase
+ xsmall<- NULL
+ xsmall$cconditional <- x$cconditional
+ zsmall$condvar <- z$condvar
repeat
{
z$n <- z$n+1
z$nit <- z$nit + 1
+ if (subphase == 1 && z$nit == 2)
+ z$time1 <- proc.time()[[3]]
if (subphase == 1 && z$nit == 11)
{
time1 <- proc.time()[[3]] - z$time1
@@ -185,7 +228,7 @@
# if (is.batch())
# {
# z$writefreq <- z$writefreq * 2 ##compensation for it
- # ## running faster with no tcl/tk
+ # ## running faster with no tcl/tk
# }
z$writefreq <- roundfreq(z$writefreq)
}
@@ -216,13 +259,27 @@
}
}
}
- zz <- x$FRAN(z, x, INIT = FALSE,...)
- if (!zz$OK)
+ if (z$int == 1)
{
- z$OK <- zz$OK
- break
+ zz <- x$FRAN(zsmall, xsmall)
+ fra <- colSums(zz$fra) - z$targets
+ if (!zz$OK)
+ {
+ z$OK <- zz$OK
+ break
+ }
}
- fra <- colSums(zz$fra) - z$targets
+ else
+ {
+ zz <- clusterCall(z$cl, usesim, zsmall, xsmall)
+ fra <- rowMeans(sapply(zz, function(x) colSums(x$fra)- z$targets))
+ zz$OK <- sapply(zz, function(x) x$OK)
+ if (!all(zz$OK))
+ {
+ z$OK <- FALSE
+ break
+ }
+ }
if (z$nit %% 2 == 1)
{
prev.fra <- fra
@@ -234,9 +291,9 @@
ac <- ifelse (z$prod0 > 1e-12, z$prod1 / z$prod0, -2)
z$maxacor <- max(-99, ac[!z$fixed]) ##note -2 > -99
z$minacor <- min(1, ac[(!z$fixed) & ac > -1.0])
+ z$ac <- ac
if (z$nit %% z$writefreq == 0)
{
- z$ac <- ac
DisplayThetaAutocor(z)
}
}
@@ -267,10 +324,11 @@
## check positivity restriction
z$positivized[z$nit, ] <- z$posj & (fchange > z$theta)
fchange <- ifelse(z$posj & (fchange > z$theta), z$theta * 0.5, fchange)
- z$theta <- z$theta - fchange
- z$thav <- z$thav + z$theta
+ zsmall$theta <- zsmall$theta - fchange
+ z$theta <- zsmall$theta
+ z$thav <- z$thav + zsmall$theta
##check for user interrupt
- ## browser()
+ ## browser()
CheckBreaks()
if (UserInterruptFlag() || UserRestartFlag() || EarlyEndPhase2Flag())
{
Modified: pkg/RSiena/R/phase3.r
===================================================================
--- pkg/RSiena/R/phase3.r 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/phase3.r 2009-07-30 21:15:49 UTC (rev 5)
@@ -11,60 +11,91 @@
# * covariance matrix
# *****************************************************************************/
##args: x model object, z control object
-phase3<- function(z, x, useCluster, noClusters, initC, ...)
+phase3<- function(z, x, ...)
{
## initialize phase 3
f <- FRANstore()
DisplayTheta(z)
z$Phase <- 3
- int <- 1
- if (useCluster)
- {
- cl <- makeCluster(rep("localhost", noClusters), type = "SOCK",
- outfile = 'cluster.out')
- clusterSetupRNG(cl, seed = rep(1, 6))
- clusterCall(cl, library, 'RSiena', character.only = TRUE)
- ## clusterExport(cl, "dllpath")
- int <- noClusters
- if (initC)
- {
- ans <- clusterCall(cl, x$FRAN, z, x, INIT=FALSE, initC = initC)
- }
- }
+ int <- z$int
+
if (x$checktime) z$ctime <- proc.time()[3]
z$sf <- matrix(0, nrow = x$n3, ncol = z$pp)
z$sf2 <- array(0, dim = c(x$n3, f$observations - 1, z$pp))
z$ssc <- array(0, dim = c(x$n3, f$observations - 1, z$pp))
z$sdf <- array(0, dim = c(x$n3, z$pp, z$pp))
- z$Deriv<-!x$FinDiff.method ## revert to original requested method for phase 3
+ ## revert to original requested method for phase 3
+ z$Deriv <- !x$FinDiff.method
if (x$FinDiff.method)
Report('Estimation of derivatives by the finite difference method.\n\n',outf)
else
- Report('Estimation of derivatives by the LR method (type 1).\n\n',outf)
+ Report('Estimation of derivatives by the LR method (type 1).\n\n', outf)
zsmall <- NULL
zsmall$theta <- z$theta
zsmall$Deriv <- z$Deriv
- zsmall$Phase<- z$Phase
+ zsmall$Phase <- z$Phase
zsmall$nit <- z$nit
xsmall<- NULL
xsmall$cconditional <- x$cconditional
zsmall$condvar <- z$condvar
- if (!x$maxlike)
+ if (!x$maxlike && !is.null(z$writefreq))
{
if (x$FinDiff.method)
z$writefreq <- z$writefreq %/% z$pp
else
z$writefreq <- z$writefreq %/% 2
- z$writefreq <-roundfreq(z$writefreq)
+ z$writefreq <- roundfreq(z$writefreq)
}
- z<-AnnouncePhase(z, x)
+ z <- AnnouncePhase(z, x)
Report('Simulated values, phase 3.\n', cf)
- for (nit in seq(1, x$n3, int))
+ if (x$n3 %% int > 0)
{
+ endNit <- x$n3 + int - x$n3 %%int
+ }
+ else
+ {
+ endNit <- x$n3
+ }
+ nits <- seq(1, endNit, int)
+ nits11 <- min(c(endNit, nits[nits >= 11]))
+ writefreq <- z$writefreq
+ if (is.null(z$writefreq))
+ {
+ z$writefreq <- 10
+ }
+ for (nit in nits)
+ {
z$nit <- nit
- if (nit <= 5 || nit==10 || nit %% z$writefreq == 0 ||
- (int > 1 && nit %% z$writefreq == 1))
+ if (is.null(writefreq))
{
+ if (nit == nits[2])
+ {
+ z$time1 <- proc.time()[[3]]
+ }
+ else if (nit == nits11)
+ {
+ time1 <- proc.time()[[3]] - z$time1
+ if (time1 > 1e-5)
+ {
+ z$writefreq <- max(1, round(20.0 / time1))
+ }
+ else
+ {
+ z$writefreq <- 20
+ }
+ if (is.batch())
+ {
+ z$writefreq <- z$writefreq * 2 ##compensation for it
+ ## running faster with no tcl/tk
+ }
+ z$writefreq <- roundfreq(z$writefreq)
+ writefreq <- z$writefreq
+ }
+ }
+ if (nit <= 5 || nit == 10 || nit %% z$writefreq == 0 ||
+ (int > 1 && nit %in% nits[seq(z$writefreq + 1, x$n3 %/% int,
+ z$writefreq)]))
+ {
if (!is.batch())
{
DisplayIteration(z)
@@ -78,53 +109,57 @@
if (nit %% z$writefreq == 0 || (int > 1 &&
nit %% z$writefreq == 1) )
{
- increment <- ifelse(nit<=5, 1, ifelse(nit <=10, 5, z$writefreq))
+ increment <- ifelse(nit <= 5, 1,
+ ifelse(nit <= 10, 5, z$writefreq))
val<- getProgressBar(z$pb)
if (x$FinDiff.method)
- val<-val+increment*(z$pp+1)
+ val <- val + increment * (z$pp + 1)
else
- val <- val+increment
+ val <- val + increment
## sort out progress bar
- z$pb <- setProgressBar(z$pb,val)
+ z$pb <- setProgressBar(z$pb, val)
if (is.batch())
- Report(c('Phase ',z$Phase,' Iteration ',nit,
- ' Progress ',round(val/z$pb$pbmax*100),'%\n'),sep='')
+ Report(c("Phase ", z$Phase, " Iteration ", nit,
+ " Progress ",
+ round(val / z$pb$pbmax * 100), "%\n"), sep='')
}
}
###############################
- z<-doPhase3it(z,x,nit,cl=cl,int=int,zsmall=zsmall,xsmall=xsmall,...)
+ z <- doPhase3it(z, x, nit, cl=z$cl, int=int, zsmall=zsmall,
+ xsmall=xsmall, ...)
##############################
- # browser()
+ ## browser()
if (!is.batch())
{
# if (nit<=3 || nit %%z$writefreq ==0)
# DisplayTheta(z)
if (nit < 10)
- Report(c(' ',nit,' ',format(z$sf[nit,],width=10,
- digits=4,nsmall=4),'\n'),cf)
+ Report(c(" ", nit, " ", format(z$sf[nit,], width=10,
+ digits=4, nsmall=4), "\n"), cf)
if (nit >= 10)
{
CheckBreaks()
## if (nit==10) set up stopkey hint Early termination of estimation
if( UserInterruptFlag())
{
- Report(c('The user asked for an early stop of the algorithm ',
- 'during phase 3.\n'),outf)
- z$Phase3Interrupt<- TRUE
- if (nit< 500)
+ Report(c("The user asked for an early stop of the algorithm ",
+ "during phase 3.\n"), outf)
+ z$Phase3Interrupt <- TRUE
+ if (nit < 500)
{
if (EarlyEndPhase2Flag())
- Report('This implies that ',outf)
+ Report('This implies that ', outf)
else
Report('This implies that the estimates are as usual,\nbut ',
outf)
Report(c('the diagnostic checks, covariance matrices and',
't-values \nare less reliable, because they are now',
- 'based on only',nit,'phase-3 iterations.\n'),outf)
+ 'based on only', nit,
+ 'phase-3 iterations.\n'), outf)
}
- z$sf<- z$sf[1:nit,,drop=FALSE]
- z$ssc<- z$ssc[1:nit,,,drop=FALSE]
- z$sdf<-z$sdf[1:nit,,,drop=FALSE]
+ z$sf <- z$sf[1:nit,,drop=FALSE]
+ z$ssc <- z$ssc[1:nit,,,drop=FALSE]
+ z$sdf <-z$sdf[1:nit,,,drop=FALSE]
break
}
if (UserRestartFlag())
@@ -132,22 +167,17 @@
}
}
}
- if (!z$OK||UserRestartFlag())
+ if (!z$OK || UserRestartFlag())
{
- if (useCluster)
- stopCluster(cl)
return(z)
}
- z$Phase3nits<- nit
+ z$Phase3nits <- nit
z<- phase3.2(z,x)
- if (useCluster)
- stopCluster(cl)
z
}
-doPhase3it<- function(z,x,nit,cl,int,zsmall,xsmall,...)
+doPhase3it<- function(z, x, nit, cl, int, zsmall, xsmall,...)
{
- z$n<- z$n + 1
if (int == 1)
{
zz <- x$FRAN(zsmall, xsmall)
@@ -157,11 +187,15 @@
z$zz <- zz
return(z)
}
+ z$n<- z$n + 1
}
else
{
- zz <- clusterCall(cl, x$FRAN, zsmall, xsmall)
- }
+ ##zz <- clusterCall(cl, simstats0c, zsmall, xsmall)
+ zz <- clusterCall(cl, usesim, zsmall, xsmall)
+ z$n <- z$n + z$int
+ # browser()
+ }
if (int == 1)
{
fra <- colSums(zz$fra)
@@ -195,9 +229,11 @@
z$sdf[nit:(nit + (int - 1)), , ] <- z$sdf0
}
else if (x$maxlike) ## as far as I can see
+ {
z$sdf[nit, , ] <- zz$dff
+ }
else
- {
+ {
if (int==1)
{
if (!is.null(zz[['sc']]))
Modified: pkg/RSiena/R/robmon.r
===================================================================
--- pkg/RSiena/R/robmon.r 2009-07-10 11:00:47 UTC (rev 4)
+++ pkg/RSiena/R/robmon.r 2009-07-30 21:15:49 UTC (rev 5)
@@ -32,13 +32,36 @@
##
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/rsiena -r 5
More information about the Rsiena-commits
mailing list