[Rsiena-commits] r242 - in pkg: RSiena RSiena/R RSiena/data RSiena/inst/doc RSiena/man RSiena/src/model/effects RSiena/tests RSienaTest RSienaTest/R RSienaTest/data RSienaTest/doc RSienaTest/inst/doc RSienaTest/man RSienaTest/src/model/effects RSienaTest/tests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Sep 10 17:55:46 CEST 2013
Author: tomsnijders
Date: 2013-09-10 17:55:45 +0200 (Tue, 10 Sep 2013)
New Revision: 242
Added:
pkg/RSienaTest/src/model/effects/AntiIsolateEffect.cpp
pkg/RSienaTest/src/model/effects/AntiIsolateEffect.h
Modified:
pkg/RSiena/DESCRIPTION
pkg/RSiena/R/phase1.r
pkg/RSiena/R/phase3.r
pkg/RSiena/R/print07Report.r
pkg/RSiena/R/printDataReport.r
pkg/RSiena/R/printInitialDescription.r
pkg/RSiena/R/siena07.r
pkg/RSiena/R/sienaeffects.r
pkg/RSiena/R/sienaprint.r
pkg/RSiena/changeLog
pkg/RSiena/data/allEffects.csv
pkg/RSiena/inst/doc/RSiena.bib
pkg/RSiena/inst/doc/RSiena_Manual.pdf
pkg/RSiena/inst/doc/RSiena_Manual.tex
pkg/RSiena/inst/doc/effects.pdf
pkg/RSiena/man/RSiena-package.Rd
pkg/RSiena/man/coCovar.Rd
pkg/RSiena/man/edit.sienaEffects.Rd
pkg/RSiena/man/effectsDocumentation.Rd
pkg/RSiena/man/getEffects.Rd
pkg/RSiena/man/hn3401.Rd
pkg/RSiena/man/iwlsm.Rd
pkg/RSiena/man/maxlikefn.Rd
pkg/RSiena/man/n3401.Rd
pkg/RSiena/man/print.sienaEffects.Rd
pkg/RSiena/man/print.sienaMeta.Rd
pkg/RSiena/man/s50.Rd
pkg/RSiena/man/s501.Rd
pkg/RSiena/man/s502.Rd
pkg/RSiena/man/s503.Rd
pkg/RSiena/man/s50a.Rd
pkg/RSiena/man/setEffect.Rd
pkg/RSiena/man/siena07.Rd
pkg/RSiena/man/siena08.Rd
pkg/RSiena/man/sienaAlgorithmCreate.Rd
pkg/RSiena/man/sienaCompositionChange.Rd
pkg/RSiena/man/sienaDataCreate.Rd
pkg/RSiena/man/sienaFit.Rd
pkg/RSiena/man/sienaGOF-auxiliary.Rd
pkg/RSiena/man/sienaGOF.Rd
pkg/RSiena/man/sienaNodeSet.Rd
pkg/RSiena/man/simstats0c.Rd
pkg/RSiena/man/tmp3.Rd
pkg/RSiena/man/tmp4.Rd
pkg/RSiena/man/varCovar.Rd
pkg/RSiena/man/varDyadCovar.Rd
pkg/RSiena/man/xtable.Rd
pkg/RSiena/src/model/effects/AllEffects.h
pkg/RSiena/src/model/effects/EffectFactory.cpp
pkg/RSiena/src/model/effects/IsolatePopEffect.cpp
pkg/RSiena/src/model/effects/IsolatePopEffect.h
pkg/RSiena/tests/parallel.Rout.save
pkg/RSienaTest/DESCRIPTION
pkg/RSienaTest/R/phase1.r
pkg/RSienaTest/R/phase3.r
pkg/RSienaTest/R/print07Report.r
pkg/RSienaTest/R/printDataReport.r
pkg/RSienaTest/R/printInitialDescription.r
pkg/RSienaTest/R/siena07.r
pkg/RSienaTest/R/sienaeffects.r
pkg/RSienaTest/R/sienaprint.r
pkg/RSienaTest/changeLog
pkg/RSienaTest/data/allEffects.csv
pkg/RSienaTest/doc/Siena_algorithms4.tex
pkg/RSienaTest/inst/doc/RSiena.bib
pkg/RSienaTest/inst/doc/RSiena_Manual.pdf
pkg/RSienaTest/inst/doc/RSiena_Manual.tex
pkg/RSienaTest/inst/doc/effects.pdf
pkg/RSienaTest/man/RSiena-package.Rd
pkg/RSienaTest/man/coCovar.Rd
pkg/RSienaTest/man/edit.sienaEffects.Rd
pkg/RSienaTest/man/effectsDocumentation.Rd
pkg/RSienaTest/man/getEffects.Rd
pkg/RSienaTest/man/hn3401.Rd
pkg/RSienaTest/man/iwlsm.Rd
pkg/RSienaTest/man/maxlikefn.Rd
pkg/RSienaTest/man/n3401.Rd
pkg/RSienaTest/man/print.sienaEffects.Rd
pkg/RSienaTest/man/print.sienaMeta.Rd
pkg/RSienaTest/man/s50.Rd
pkg/RSienaTest/man/s501.Rd
pkg/RSienaTest/man/s502.Rd
pkg/RSienaTest/man/s503.Rd
pkg/RSienaTest/man/s50a.Rd
pkg/RSienaTest/man/setEffect.Rd
pkg/RSienaTest/man/siena07.Rd
pkg/RSienaTest/man/siena08.Rd
pkg/RSienaTest/man/sienaAlgorithmCreate.Rd
pkg/RSienaTest/man/sienaBayes.Rd
pkg/RSienaTest/man/sienaCompositionChange.Rd
pkg/RSienaTest/man/sienaDataCreate.Rd
pkg/RSienaTest/man/sienaFit.Rd
pkg/RSienaTest/man/sienaGOF-auxiliary.Rd
pkg/RSienaTest/man/sienaGOF.Rd
pkg/RSienaTest/man/sienaNodeSet.Rd
pkg/RSienaTest/man/simstats0c.Rd
pkg/RSienaTest/man/tmp3.Rd
pkg/RSienaTest/man/tmp4.Rd
pkg/RSienaTest/man/varCovar.Rd
pkg/RSienaTest/man/varDyadCovar.Rd
pkg/RSienaTest/man/xtable.Rd
pkg/RSienaTest/src/model/effects/AllEffects.h
pkg/RSienaTest/src/model/effects/EffectFactory.cpp
pkg/RSienaTest/tests/parallel.Rout.save
Log:
Correction of an error in the dolby Option. Further, three new effects and some minor changes.
Modified: pkg/RSiena/DESCRIPTION
===================================================================
--- pkg/RSiena/DESCRIPTION 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/DESCRIPTION 2013-09-10 15:55:45 UTC (rev 242)
@@ -1,8 +1,8 @@
Package: RSiena
Type: Package
Title: Siena - Simulation Investigation for Empirical Network Analysis
-Version: 1.1-241
-Date: 2013-08-23
+Version: 1.1-242
+Date: 2013-09-10
Author: Ruth Ripley, Krists Boitmanis, Tom A.B. Snijders
Depends: R (>= 2.15.0)
Imports: Matrix
Modified: pkg/RSiena/R/phase1.r
===================================================================
--- pkg/RSiena/R/phase1.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/phase1.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -145,7 +145,7 @@
z$regrCor <- rep(0, z$pp)
if (!is.null(z$ssc))
{
- scores <- apply(z$ssc, c(1,3), mean)
+ scores <- apply(z$ssc, c(1,3), sum)
for (i in 1:z$pp)
{
if (var(scores[,i]) > 0)
@@ -217,7 +217,7 @@
fchange <- as.vector(dinv %*% z$mnfra)
z$dinv <- dinv
}
- if (!x$diagg)
+ if (!x$diagg)
{
# Partial diagonalization of derivative matrix
# for use if 0 < x$diagonalize < 1.
Modified: pkg/RSiena/R/phase3.r
===================================================================
--- pkg/RSiena/R/phase3.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/phase3.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -103,114 +103,121 @@
'with better initial values or a reduced model.\n',
'(Check that you entered the data properly!)\n'), outf)
}
- Heading(2, outf, c('End of stochastic approximation algorithm, phase ',
- z$Phase, '.'))
- Report(c('Total of', z$n,'iterations.\n'), outf)
- Report(c('Parameter estimates based on', z$n - z$Phase3nits,
- 'iterations,\n'), outf)
- if (!x$maxlike && z$cconditional)
+ if ((x$nsub == 0)&(x$simOnly))
{
- Report(c('basic rate parameter',
- c('', 's')[as.integer(z$f$observations > 2) + 1],
- ' as well as \n'), sep='', outf)
+ # nothing
}
- Report(c('convergence diagnostics, covariance and derivative matrices ',
- 'based on ', z$Phase3nits, ' iterations.\n\n'), sep='', outf)
- Report('Information for convergence diagnosis.\n', outf)
- Report(c('Averages, standard deviations, ',
- 'and t-ratios for deviations from targets:\n'), sep='', outf)
- # Report(c(date(),'\n'),bof)
- if (x$maxlike)
+ else
{
- Report('\nMaximum Likelihood estimation.', bof)
- }
- else
- {
- if (z$cconditional)
- {
- Report('\nconditional moment estimation.', bof)
+ Heading(2, outf, c('End of stochastic approximation algorithm, phase ',
+ z$Phase, '.'))
+ Report(c('Total of', z$n,'iterations.\n'), outf)
+ Report(c('Parameter estimates based on', z$n - z$Phase3nits,
+ 'iterations,\n'), outf)
+ if (!x$maxlike && z$cconditional)
+ {
+ Report(c('basic rate parameter',
+ c('', 's')[as.integer(z$f$observations > 2) + 1],
+ ' as well as \n'), sep='', outf)
+ }
+ Report(c('convergence diagnostics, covariance and derivative matrices ',
+ 'based on ', z$Phase3nits, ' iterations.\n\n'), sep='', outf)
+ Report('Information for convergence diagnosis.\n', outf)
+ Report(c('Averages, standard deviations, ',
+ 'and t-ratios for deviations from targets:\n'), sep='', outf)
+ # # Report(c(date(),'\n'),bof)
+ if (x$maxlike)
+ {
+ Report('\nMaximum Likelihood estimation.', bof)
+ }
+ else
+ {
+ if (z$cconditional)
+ {
+ Report('\nconditional moment estimation.', bof)
- }
- else
- {
- Report('\nunconditional moment estimation.', bof)
- }
+ }
+ else
+ {
+ Report('\nunconditional moment estimation.', bof)
+ }
+ }
+ Report('\nInformation for convergence diagnosis.\n', bof)
+ Report(c('Averages, standard deviations, ',
+ 'and t-ratios for deviations from targets:\n'), bof, sep='')
+ ##calculate t-ratios
+ dmsf <- diag(z$msf)
+ sf <- colMeans(z$sf)
+ # TS: I wonder why "use" and "use2" are the names in the following lines;
+ # the coordinates with "use" are <<not>> used.
+ use <- dmsf < 1e-20 * z$scale * z$scale
+ use2 <- abs(sf) < 1e-10 * z$scale
+ dmsf[use] <- 1e-20 * z$scale[use] * z$scale[use]
+ tstat <- rep(NA, z$pp)
+ tstat[!use]<- sf[!use] / sqrt(dmsf[!use])
+ tstat[use & use2] <- 0
+ tstat[use & !use2] <- 999
+ z$tstat <- tstat
+ # tconv.max = Maximum value of t-ratio for convergence,
+ # for any linear combination.
+ z$tconv.max <- NA
+ if (sum(!z$fixed) > 0)
+ {
+ mean.dev <- colSums(z$sf)[!z$fixed]/dim(z$sf)[1]
+ cov.dev <- z$msf[!z$fixed,!z$fixed]
+ if (inherits(try(thisproduct <- solve(cov.dev, mean.dev), silent=TRUE),
+ "try-error"))
+ {
+ Report('Maximum t-ratio for convergence not computable.\n', outf)
+ }
+ else
+ {
+ z$tconv.max <- sqrt(t(mean.dev) %*% thisproduct)
+ }
+ }
+ mymess1 <- paste(format(1:z$pp,width=3), '. ',
+ format(round(sf, 4), width=8, nsmall=4), ' ',
+ format(round(sqrt(dmsf), 4) ,width=8, nsmall=4), ' ',
+ format(round(tstat, 4), width=8, nsmall=3), sep='')
+ mymess2 <- c('', ' (fixed parameter)')[as.numeric(z$fixed) + 1]
+ mymess <- paste(mymess1, mymess2)
+ PrtOutMat(as.matrix(mymess), outf)
+ PrtOutMat(as.matrix(mymess1), bof)
+ ## Report(mymess1, bof, fill=80)
+ tmax <- max(abs(tstat)[!z$fixed & !z$BasicRateFunction & z$resist > 0.9])
+ z$tconv <- tstat
+ error <- (abs(tmax) > 4.0 / sqrt(z$Phase3nits)) && (abs(tmax) > 0.3)
+ if (tmax >= 0.4 & !z$error)
+ {
+ z$error <- TRUE
+ }
+ Report('Good convergence is indicated by the t-ratios ', outf)
+ if (any(z$fixed))
+ {
+ Report('of non-fixed parameters ', outf)
+ }
+ Report('being close to zero.\n', outf)
+ if (z$Phase3nits < 100)
+ {
+ Report(c('(Since the diagnostic checks now are based only on ',
+ z$Phase3nits,
+ ' iterations,', '\nthey are not reliable.)\n'), sep='', outf)
+ }
+ if (error) ## also test subphase here but not relevant to phase 3, I think
+ {
+ Report('One or more of the t-statistics are rather large.\n', outf)
+ if (tmax > 0.5)
+ {
+ Report('Convergence of the algorithm is doubtful.\n', outf)
+ }
+ ## removed repfortotal loop possibility here as not functioning now
+ if (z$Phase3nits <= 50)
+ {
+ Report(c('Note that the standard deviations are based on',
+ 'few simulations.\n'), outf)
+ }
+ }
}
- Report('\nInformation for convergence diagnosis.\n', bof)
- Report(c('Averages, standard deviations, ',
- 'and t-ratios for deviations from targets:\n'), bof, sep='')
- ##calculate t-ratios
- dmsf <- diag(z$msf)
- sf <- colMeans(z$sf)
- # TS: I wonder why "use" and "use2" are the names in the following lines;
- # the coordinates with "use" are <<not>> used.
- use <- dmsf < 1e-20 * z$scale * z$scale
- use2 <- abs(sf) < 1e-10 * z$scale
- dmsf[use] <- 1e-20 * z$scale[use] * z$scale[use]
- tstat <- rep(NA, z$pp)
- tstat[!use]<- sf[!use] / sqrt(dmsf[!use])
- tstat[use & use2] <- 0
- tstat[use & !use2] <- 999
- z$tstat <- tstat
- # tconv.max = Maximum value of t-ratio for convergence,
- # for any linear combination.
- z$tconv.max <- NA
- if (sum(!z$fixed) > 0)
- {
- mean.dev <- colSums(z$sf)[!z$fixed]/dim(z$sf)[1]
- cov.dev <- z$msf[!z$fixed,!z$fixed]
- if (inherits(try(thisproduct <- solve(cov.dev, mean.dev), silent=TRUE),
- "try-error"))
- {
- Report('Maximum t-ratio for convergence not computable.\n', outf)
- }
- else
- {
- z$tconv.max <- sqrt(t(mean.dev) %*% thisproduct)
- }
- }
- mymess1 <- paste(format(1:z$pp,width=3), '. ',
- format(round(sf, 4), width=8, nsmall=4), ' ',
- format(round(sqrt(dmsf), 4) ,width=8, nsmall=4), ' ',
- format(round(tstat, 4), width=8, nsmall=3), sep='')
- mymess2 <- c('', ' (fixed parameter)')[as.numeric(z$fixed) + 1]
- mymess <- paste(mymess1, mymess2)
- PrtOutMat(as.matrix(mymess), outf)
- PrtOutMat(as.matrix(mymess1), bof)
- ## Report(mymess1, bof, fill=80)
- tmax <- max(abs(tstat)[!z$fixed & !z$BasicRateFunction & z$resist > 0.9])
- z$tconv <- tstat
- error <- (abs(tmax) > 4.0 / sqrt(z$Phase3nits)) && (abs(tmax) > 0.3)
- if (tmax >= 0.4 & !z$error)
- {
- z$error <- TRUE
- }
- Report('Good convergence is indicated by the t-ratios ', outf)
- if (any(z$fixed))
- {
- Report('of non-fixed parameters ', outf)
- }
- Report('being close to zero.\n', outf)
- if (z$Phase3nits < 100)
- {
- Report(c('(Since the diagnostic checks now are based only on ',
- z$Phase3nits,
- ' iterations,', '\nthey are not reliable.)\n'), sep='', outf)
- }
- if (error) ## also test subphase here but not relevant to phase 3, I think
- {
- Report('One or more of the t-statistics are rather large.\n', outf)
- if (tmax > 0.5)
- {
- Report('Convergence of the algorithm is doubtful.\n', outf)
- }
- ## removed repfortotal loop possibility here as not functioning now
- if (z$Phase3nits <= 50)
- {
- Report(c('Note that the standard deviations are based on',
- 'few simulations.\n'), outf)
- }
- }
if (x$maxlike)
{
Report('Autocorrelations during phase 3 : \n', outf)
@@ -260,7 +267,7 @@
if (inherits(try(cov <- solve(dfrac), silent=TRUE),"try-error"))
{
Report('Noninvertible estimated covariance matrix : \n', outf)
- errorMessage.cov <- 'Noninvertible estimated covariance matrix'
+ errorMessage.cov <- '***Warning: Noninvertible estimated covariance matrix ***'
cov <- NULL
}
}
@@ -271,7 +278,7 @@
error <- FALSE
if (inherits(try(msfinv <- solve(z$msfc), silent=TRUE), "try-error"))
{
- Report('Covariance matrix not positive definite: \n', outf)
+ Report('*** Warning: Covariance matrix not positive definite *** \n', outf)
if (any(z$fixed || any(z$newfixed)))
{
Report(c('(This may be unimportant, and related to the fact\n',
@@ -281,7 +288,7 @@
{
Report(c('This may mean that the reported standard errors ',
'are invalid.\n'), outf)
- errorMessage.cov <- 'Noninvertible estimated covariance matrix'
+ errorMessage.cov <- '*** Warning: Noninvertible estimated covariance matrix ***'
}
z$msfinv <- NULL
}
@@ -310,6 +317,7 @@
CalculateDerivative3<- function(z,x)
{
z$mnfra <- colMeans(z$sf)
+ estMeans <- z$mnfra + z$targets
if (z$FinDiff.method || x$maxlike)
{
dfra <- t(as.matrix(Reduce("+", z$sdf) / length(z$sdf)))
@@ -328,18 +336,24 @@
Report(c("Warning: diagonal element(s)", sub,
" of derivative matrix < 0\n"), cf)
}
- scores <- apply(z$ssc, c(1,3), mean)
+ scores <- apply(z$ssc, c(1,3), sum) # z$nit by z$pp matrix
for (i in 1:z$pp)
{
- if (var(scores[,i]) > 0)
+ if ((var(scores[,i]) > 0)&(var(z$sf[,i]) > 0))
{
z$regrCoef[i] <- cov(z$sf[,i], scores[,i])/var(scores[,i])
z$regrCor[i] <- cor(z$sf[,i], scores[,i])
}
+ }
+ if (x$dolby)
+ {
+ mean.scores <- colMeans(scores)
+ estMeans <- estMeans - (z$regrCoef * colMeans(scores))
}
Report('Correlations between scores and statistics:\n', cf)
PrtOutMat(format(as.matrix(t(z$regrCor)), digits = 2, nsmall = 2), cf)
}
+ z$estMeans <- estMeans
z$diver <- rep(FALSE, z$pp)
if (z$AllUserFixed & any(abs(diag(dfra)) < 1e-6))
{
@@ -379,7 +393,7 @@
{
Report(c('Warning. After phase 3, derivative matrix non-invertible',
'even with a ridge.\n'), cf)
- z$errorMessage.cov <- 'Noninvertible derivative matrix'
+ z$errorMessage.cov <- '*** Warning: Noninvertible derivative matrix ***'
fchange <- 0
z$dinv <- z$dfrac
z$dinv[,] <- 999
Modified: pkg/RSiena/R/print07Report.r
===================================================================
--- pkg/RSiena/R/print07Report.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/print07Report.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -11,197 +11,230 @@
##@PrintReport siena07 Print report
PrintReport <- function(z, x)
{
- types <- attr(z$f, "types")
- Report('\n\n', outf)
+ types <- attr(z$f, "types")
+ Report('\n\n', outf)
if ((x$nsub == 0)&(x$simOnly))
{
Heading(2, outf, "Simulation Results.")
}
else
- {
+ {
Heading(2, outf, "Estimation Results.")
}
if (!z$OK)
- {
- Report("Error end of estimation algorithm", outf)
- }
- else
- {
- Report("Regular end of estimation algorithm.\n", outf)
- Report(c("Total of", z$n, "iteration steps.\n\n"), outf)
- Report(c("Total of", z$n, "iteration steps.\n\n"), bof)
- if (x$simOnly)
- {
- Heading(3, outf, "Parameter values")
- }
- else
+ {
+ Report("Error end of estimation algorithm", outf)
+ }
+ else
+ {
+ if (x$simOnly)
+ {
+ Report("Regular end of simulation algorithm.\n", outf)
+ Report(c("Total of", z$n, "iteration steps.\n\n"), outf)
+ Report(c("Total of", z$n, "iteration steps.\n\n"), bof)
+ Heading(3, outf, "Parameter values")
+ }
+ else
+ {
+ Report("Regular end of estimation algorithm.\n", outf)
+ Report(c("Total of", z$n, "iteration steps.\n\n"), outf)
+ Report(c("Total of", z$n, "iteration steps.\n\n"), bof)
+ Heading(3, outf, "Estimates and standard errors")
+ }
+ Heading(3, bof, "Estimates and standard errors")
+ if (z$cconditional) ## deal with rate parameter
+ {
+ Report('Rate parameters:\n', outf)
+ Report('Rate parameters:\n', bof)
+ if (length(z$rate) == 1)
{
- Heading(3, outf, "Estimates and standard errors")
+ if (length(attr(z$f,'netnames')) == 1)
+ {
+ Report(format(' 0. Rate parameter', width = 43), outf)
+ Report(format(' 0. Rate parameter', width = 43), bof)
+ }
+ else
+ {
+ Report(format(' 0. Rate parameter of conditioning variable',
+ width = 43), outf)
+ Report(format(' 0. Rate parameter of conditioning variable',
+ width = 43), bof)
+ }
+ Report(format(round(z$rate[1], digits = 4), width = 9), outf)
+ Report(format(round(z$rate[1], digits = 4), width = 9), bof)
+ Report(c(' (', format(round(z$vrate[1], digits = 4),
+ width = 9), ')\n'), sep = '', outf)
+ Report(c(' (', format(round(z$vrate[1], digits = 4),
+ width = 9), ')\n'), sep = '', bof)
}
- Heading(3, bof, "Estimates and standard errors")
- if (z$cconditional) ## deal with rate parameter
- {
- Report('Rate parameters:\n', outf)
- Report('Rate parameters:\n', bof)
- if (length(z$rate) == 1)
- {
- if (length(attr(z$f,'netnames')) == 1)
- {
- Report(format(' 0. Rate parameter', width = 43), outf)
- Report(format(' 0. Rate parameter', width = 43), bof)
- }
- else
- {
- Report(format(' 0. Rate parameter of conditioning variable',
- width = 43), outf)
- Report(format(' 0. Rate parameter of conditioning variable',
- width = 43), bof)
- }
- Report(format(round(z$rate[1], digits = 4), width = 9), outf)
- Report(format(round(z$rate[1], digits = 4), width = 9), bof)
- Report(c(' (', format(round(z$vrate[1], digits = 4),
- width = 9), ')\n'), sep = '', outf)
- Report(c(' (', format(round(z$vrate[1], digits = 4),
- width = 9), ')\n'), sep = '', bof)
- }
- else ## observations > 2
- {
- nn <- length(z$rate)
- nnstr <- format(as.character(1:nn),width=2,justify='left')
- if (z$f$nDepvars == 1 || is.null(z$f$nDepvars))
- {
- tmp <- paste(' 0.', nnstr, ' Rate parameter period ',
- 1:nn, ' ',
- format(round(z$rate,4),width=9),
- ' (',format(round(z$vrate,4),width=9),
- ')\n', sep = '')
- } else{
- tmp <- paste(' 0.', nnstr,
- 'Rate parameter cond. variable period ',
- 1:nn, ' ',
- format(round(z$rate,4),width=9),
- ' (',format(round(z$vrate,4),width=9),
- ')\n', sep='')
- }
- Report(tmp, outf, sep='')
- Report(tmp, bof, sep='')
- Report('\nOther parameters:\n', outf)
- Report('\nOther parameters:\n', bof)
- }
- }
- nBehavs <- sum(types == "behavior")
- nNetworks <- length(types) - nBehavs
- if (nBehavs > 0 && nNetworks > 0)
- {
- Report("Network Dynamics\n", outf)
- }
- if (!x$simOnly)
+ else ## observations > 2
{
- ses <- ifelse(diag(z$covtheta) >= 0.0 | z$fixed,
- paste(' (', sprintf("%9.4f",sqrt(diag(z$covtheta))),
- ')', sep=''), ' ---')
- if (!all(z$fixed))
+ nn <- length(z$rate)
+ nnstr <- format(as.character(1:nn),width=2,justify='left')
+ if (z$f$nDepvars == 1 || is.null(z$f$nDepvars))
{
- ses[z$fixed] <- ' ( fixed )'
+ tmp <- paste(' 0.', nnstr, ' Rate parameter period ',
+ 1:nn, ' ',
+ format(round(z$rate,4),width=9),
+ ' (',format(round(z$vrate,4),width=9),
+ ')\n', sep = '')
+ } else{
+ tmp <- paste(' 0.', nnstr,
+ 'Rate parameter cond. variable period ',
+ 1:nn, ' ',
+ format(round(z$rate,4),width=9),
+ ' (',format(round(z$vrate,4),width=9),
+ ')\n',sep='')
}
- theta <- ifelse(diag(z$covtheta) >= 0.0 | z$fixed,
- sprintf("%9.4f",z$theta),
- ' ---')
- }
- else
+ Report(tmp, outf, sep='')
+ Report(tmp, bof, sep='')
+ }
+ Report('\nOther parameters:\n', outf)
+ Report('\nOther parameters:\n', bof)
+ }
+ nBehavs <- sum(types == "behavior")
+ nNetworks <- length(types) - nBehavs
+ if (nBehavs > 0 && nNetworks > 0)
+ {
+ Report("Network Dynamics\n", outf)
+ }
+ if (!x$simOnly)
+ {
+ ses <- ifelse(diag(z$covtheta) >= 0.0 | z$fixed,
+ paste(' (', sprintf("%9.4f",sqrt(diag(z$covtheta))),
+ ')', sep=''), ' ---')
+ if (!all(z$fixed))
{
- theta <- sprintf("%9.4f",z$theta)
+ ses[z$fixed] <- ' ( fixed )'
}
- if (nBehavs > 0)
+ theta <- ifelse(diag(z$covtheta) >= 0.0 | z$fixed,
+ sprintf("%9.4f",z$theta),
+ ' ---')
+ }
+ else
+ {
+ theta <- sprintf("%9.4f",z$theta)
+ }
+ if (nBehavs > 0)
+ {
+ behEffects <-
+ z$requestedEffects[z$requestedEffects$netType == 'behavior',]
+ behNames <- unique(behEffects$name)
+ if (nBehavs > 1)
{
- behEffects <-
- z$requestedEffects[z$requestedEffects$netType == 'behavior',]
- behNames <- unique(behEffects$name)
- if (nBehavs > 1)
- {
- behEffects$effectName <- paste('<',
- (1:nBehavs)[match(behEffects$name,
- behNames)],
- '> ', behEffects$effectName,
- sep='')
- z$requestedEffects$effectName[z$requestedEffects$netType==
- 'behavior'] <-
- behEffects$effectName
- }
- }
- typesp <- ifelse (z$requestedEffects$type %in% c("eval", "rate"),
- ": ", ": ")
- typetxt <- ifelse (z$requestedEffects$type == "creation", "creat",
- z$requestedEffects$type )
- tmp <- paste(typetxt, typesp, z$requestedEffects$effectName, sep = '')
- if (x$simOnly)
+ behEffects$effectName <- paste('<',
+ (1:nBehavs)[match(behEffects$name,
+ behNames)],
+ '> ', behEffects$effectName,
+ sep='')
+ z$requestedEffects$effectName[z$requestedEffects$netType==
+ 'behavior'] <-
+ behEffects$effectName
+ }
+ }
+ typesp <- ifelse (z$requestedEffects$type %in% c("eval", "rate"),
+ ": ", ": ")
+ typetxt <- ifelse (z$requestedEffects$type == "creation", "creat",
+ z$requestedEffects$type )
+ tmp <- paste(typetxt, typesp, z$requestedEffects$effectName, sep = '')
+ if (x$simOnly)
+ {
+ ses <- rep(' ', z$pp)
+ }
+ tmp <- paste(sprintf("%2d", 1:length(z$requestedEffects$effectName)),
+ '. ', format(substr(tmp, 1, 50), width=50),
+ theta, ses, '\n', sep='', collapse = '')
+ if (nBehavs > 0 && nNetworks > 0)
+ {
+ nNetworkEff <- nrow(z$requestedEffects) - nrow(behEffects)
+ tmpstr <- paste(nNetworkEff + 1, '. ', sep='')
+ tmpsub <- regexpr(tmpstr, tmp, fixed=TRUE)
+ tmp1 <- substring(tmp, 1, tmpsub - 2)
+ tmp2 <- substring(tmp, tmpsub - 1, nchar(tmp))
+ Report(tmp1, outf)
+ Report(tmp1, bof)
+ Report('\nBehavior Dynamics\n', outf)
+ Report(tmp2, outf)
+ Report(tmp2, bof)
+ }
+ else
+ {
+ Report(tmp, outf)
+ Report(tmp, bof)
+ }
+ Report('\n', outf)
+ Report('\n', bof)
+ if (z$cconditional && length(attr(z$f, 'netnames')) > 1)
+ {
+ Report(c('For conditional estimation, ',
+ 'the standard errors of rate parameters\n',
+ 'not used for conditioning are unreliable.'), outf)
+ }
+ if (x$simOnly)
+ {
+ Heading(3, outf, "Simulated statistics")
+ Report('Estimated means and standard deviations, standard errors of the mean \n', bof)
+ Report('Estimated means and standard deviations, standard errors of the mean \n', outf)
+ dmsf <- diag(z$msf)
+ sf <- colMeans(z$sf)
+ mean.stats <- colMeans(z$sf) + z$targets
+ cov.dev <- z$msf
+ sem <- sqrt(dmsf/dim(z$sf)[1])
+ if (x$dolby)
{
- ses <- rep(' ', z$pp)
- }
- tmp <- paste(sprintf("%2d", 1:length(z$requestedEffects$effectName)),
- '. ', format(substr(tmp, 1, 50), width=50),
- theta, ses, '\n', sep='', collapse = '')
- if (nBehavs > 0 && nNetworks > 0)
- {
- nNetworkEff <- nrow(z$requestedEffects) - nrow(behEffects)
- tmpstr <- paste(nNetworkEff + 1, '. ', sep='')
- tmpsub <- regexpr(tmpstr, tmp, fixed=TRUE)
- tmp1 <- substring(tmp, 1, tmpsub - 2)
- tmp2 <- substring(tmp, tmpsub - 1, nchar(tmp))
- Report(tmp1, outf)
- Report(tmp1, bof)
- Report('\nBehavior Dynamics\n', outf)
- Report(tmp2, outf)
- Report(tmp2, bof)
- }
- else
- {
- Report(tmp, outf)
- Report(tmp, bof)
- }
- Report('\n', outf)
- Report('\n', bof)
- if (z$cconditional && length(attr(z$f, 'netnames')) > 1)
- {
- Report(c('For conditional estimation, ',
- 'the standard errors of rate parameters\n',
- 'not used for conditioning are unreliable.'), outf)
- }
- if (!x$simOnly)
+ scores <- apply(z$ssc, c(1,3), sum) # z$nit by z$pp matrix
+ mean.scores <- colMeans(scores)
+ mean.stats <- mean.stats - (z$regrCoef * mean.scores)
+ sem <- sem*sqrt(1 - (z$regrCor)^2)
+ }
+ mymess1 <- paste(format(1:z$pp,width=3), '. ',
+ format(z$requestedEffects$functionName, width = 56),
+ format(round(mean.stats, 3), width=8, nsmall=3), ' ',
+ format(round(sqrt(dmsf), 3) ,width=8, nsmall=3), ' ',
+ format(round(sem, 4) ,width=8, nsmall=4), sep='')
+ PrtOutMat(as.matrix(mymess1), outf)
+ PrtOutMat(as.matrix(mymess1), bof)
+ if (x$dolby)
{
- Heading(3, outf, "Covariance matrices")
- if (any(z$fixed))
- {
- Report(c('(Values of the covariance matrix of estimates\n',
- ' are meaningless for the fixed parameters.)\n\n'),
- outf)
- }
- Report(c("Covariance matrix of estimates",
- "(correlations below diagonal):\n"), outf)
- covcor <- z$covtheta
- correl <- z$covtheta/sqrt(diag(z$covtheta))[row(z$covtheta)]/
- sqrt(diag(z$covtheta))[col(z$covtheta)]
- covcor[lower.tri(covcor)] <- correl[lower.tri(correl)]
- PrtOutMat(format(round(covcor, digits = 3), width = 10), outf)
- Report(c('Derivative matrix of expected statistics X by',
- 'parameters and\n'), outf)
- Report(c("covariance/correlation matrix of X can be found using\n",
- "summary(ans) within R,",
- " or by using the 'verbose' option in Siena07.\n "),
- sep = "", outf)
- Report(c("Derivative matrix of expected statistics X by",
- "parameters:\n\n"), lf)
- PrtOutMat(z$dfrac, lf)
- Report('Covariance matrix of X (correlations below the diagonal):\n',
- lf)
- covcor <- z$msf
- correl <- z$msf/sqrt(diag(z$msf))[row(z$msf)]/
- sqrt(diag(z$msf))[col(z$msf)]
- covcor[lower.tri(covcor)] <- correl[lower.tri(correl)]
- PrtOutMat(format(round(covcor, digits = 3), width = 10), lf)
- Report('\n', outf)
- Report('\n', lf)
- }
+ Report('Standard errors of the mean are less than s.d./n \n', outf)
+ Report('because of regression on scores (Dolby option). \n', outf)
+ }
}
+ else
+ {
+ Heading(3, outf, "Covariance matrices")
+ if (any(z$fixed))
+ {
+ Report(c('(Values of the covariance matrix of estimates\n',
+ ' are meaningless for the fixed parameters.)\n\n'),
+ outf)
+ }
+ Report(c("Covariance matrix of estimates",
+ "(correlations below diagonal):\n"), outf)
+ covcor <- z$covtheta
+ correl <- z$covtheta/sqrt(diag(z$covtheta))[row(z$covtheta)]/
+ sqrt(diag(z$covtheta))[col(z$covtheta)]
+ covcor[lower.tri(covcor)] <- correl[lower.tri(correl)]
+ PrtOutMat(format(round(covcor, digits = 3), width = 10), outf)
+ Report(c('Derivative matrix of expected statistics X by',
+ 'parameters and\n'), outf)
+ Report(c("covariance/correlation matrix of X can be found using\n",
+ "summary(ans) within R,",
+ " or by using the 'verbose' option in Siena07.\n "),
+ sep = "", outf)
+ Report(c("Derivative matrix of expected statistics X by",
+ "parameters:\n\n"), lf)
+ PrtOutMat(z$dfrac, lf)
+ Report('Covariance matrix of X (correlations below the diagonal):\n',
+ lf)
+ covcor <- z$msf
+ correl <- z$msf/sqrt(diag(z$msf))[row(z$msf)]/
+ sqrt(diag(z$msf))[col(z$msf)]
+ covcor[lower.tri(covcor)] <- correl[lower.tri(correl)]
+ PrtOutMat(format(round(covcor, digits = 3), width = 10), lf)
+ Report('\n', outf)
+ Report('\n', lf)
+ }
+ }
}
Modified: pkg/RSiena/R/printDataReport.r
===================================================================
--- pkg/RSiena/R/printDataReport.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/printDataReport.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -132,8 +132,16 @@
}
if (z$FinDiff.method)
{
+ if ((x$nsub == 0)&(x$simOnly))
+ {
+ Report(c('Derivatives are estimated with the finite difference',
+ 'method.\n'), outf)
+ }
+ else
+ {
Report(c('Standard errors are estimated with the finite difference',
'method.\n'), outf)
+ }
if (!x$FinDiff.method)
{
Report("Note that the option requested has been over-ridden\n",
@@ -142,9 +150,21 @@
}
else
{
- Report(c('Standard errors are estimated with the likelihood ratio',
+ if ((x$nsub == 0)&(x$simOnly))
+ {
+ Report(c('Derivatives are estimated with the likelihood ratio',
'method.\n'), outf)
+ }
+ else
+ {
+ Report(c('Standard errors are estimated with the likelihood ratio',
+ 'method.\n'), outf)
+ }
}
+ if (x$dolby)
+ {
+ Report('Dolby method (regression on scores) is used.\n', outf)
+ }
## any max degree restrictions?
if (any(x$MaxDegree > 0))
{
@@ -181,10 +201,20 @@
}
}
}
- Report(sprintf("Initial value of gain parameter is %10.7f.\n", x$firstg),
- outf)
- Report(sprintf("Number of subphases in Phase 2 is %d.\n\n", x$nsub), outf)
- Report('Initial parameter values are \n', outf)
+ if ((x$nsub == 0)&(x$simOnly))
+ {
+ Report('Parameter values are \n', outf)
+ }
+ else
+ {
+ Report(sprintf("Initial value of gain parameter is %10.7f.\n",
+ x$firstg), outf)
+ Report(sprintf("Reduction factor for gain parameter is %10.7f.\n",
+ x$reduceg), outf)
+ Report(sprintf("Number of subphases in Phase 2 is %d.\n\n", x$nsub),
+ outf)
+ Report('Initial parameter values are \n', outf)
+ }
if (z$cconditional)
{
if (observations == 1)
Modified: pkg/RSiena/R/printInitialDescription.r
===================================================================
--- pkg/RSiena/R/printInitialDescription.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/printInitialDescription.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -232,7 +232,7 @@
}
if (!bipartite)
{
- Report("Dyad Counts:\n", outf)
+ Report("Directed dyad Counts:\n", outf)
if (valmin == 0 & valmax == 1)
{
Report(" observation total mutual asymm. null\n",
Modified: pkg/RSiena/R/siena07.r
===================================================================
--- pkg/RSiena/R/siena07.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/siena07.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -200,7 +200,14 @@
format(as.Date(packageDescription(pkgname, fields = "Date")),
"%d %b %y"), ")",
revision, "\n\n"), sep = '', outf )
+ if (z$x$simOnly)
+ {
+ Heading(1, outf, "Simulations.")
+ }
+ else
+ {
Heading(1, outf, "Estimation by stochastic approximation algorithm.")
+ }
if (is.null(seed))
{
Report("Random initialization of random number stream.\n", outf)
Modified: pkg/RSiena/R/sienaeffects.r
===================================================================
--- pkg/RSiena/R/sienaeffects.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/sienaeffects.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -12,6 +12,10 @@
type="eval", interaction1="", interaction2="",
character=FALSE)
{
+ if (!inherits(myeff, 'sienaEffects'))
+ {
+ stop("The first argument is not of class <sienaEffects>.")
+ }
if (character)
{
dots <- sapply(list(...), function(x)x)
@@ -48,8 +52,9 @@
}
else
{
- print.data.frame(myeff[use, c("name", "shortName", "type",
- "interaction1", "interaction2", "include")])
+# print.data.frame(myeff[use, c("name", "shortName", "type",
+# "interaction1", "interaction2", "include")])
+ print.sienaEffects(myeff[use,])
}
myeff
}
@@ -244,8 +249,10 @@
myeff[use, "test"] <- test
myeff[use, "initialValue"] <- initialValue
myeff[use, "timeDummy"] <- timeDummy
- print.data.frame(myeff[use, c("name", "shortName", "type", "interaction1",
- "interaction2", "include", "parm", "fix", "test",
- "initialValue", "timeDummy", "period", "group")])
+# print.data.frame(myeff[use, c("name", "shortName", "type", "interaction1",
+# "interaction2", "include", "parm", "fix", "test",
+# "initialValue", "timeDummy", "period", "group")])
+ print.sienaEffects(myeff[use,])
+
myeff
}
Modified: pkg/RSiena/R/sienaprint.r
===================================================================
--- pkg/RSiena/R/sienaprint.r 2013-08-23 15:23:15 UTC (rev 241)
+++ pkg/RSiena/R/sienaprint.r 2013-09-10 15:55:45 UTC (rev 242)
@@ -202,23 +202,47 @@
}
else
{
+ if(x$x$simOnly)
+ {
+ cat("Parameter values\n\n")
+ }
+ else
+ {
cat("Estimates, standard errors and convergence t-ratios\n\n")
+ }
tmp <- sienaFitThetaTable(x, tstat=tstat)
mydf <- tmp$mydf
mymat <- as.matrix(mydf)
mymat[, 'value'] <- format(round(mydf$value, digits=4))
- mymat[, 'se'] <- format(round(mydf$se, digits=4))
- mymat[, 'tstat'] <- paste(" ", format(round(mydf$tstat, digits=4)))
- mymat[is.na(mydf$tstat), 'tstat'] <- ' '
+ if(x$x$simOnly)
+ {
+ mymat[, 'se'] <- ' '
+ mymat[, 'tstat'] <- ' '
+ }
+ else
+ {
+ mymat[, 'se'] <- format(round(mydf$se, digits=4))
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/rsiena -r 242
More information about the Rsiena-commits
mailing list