## NULLThis website converted the following original .R scripts into .rmd files.
Please visit GitHub for the latest .R files.
Specific questions with respect to the .rmd files can be addressed to: Jochem Tolsma.
For questions on RSiena please visit the designated GitHub page.
to do: conver to .rmd
################################################################################ ----
################################################################################ Rscript04SienaBehaviour.R:
################################################################################ a script for the
################################################################################ introduction to
################################################################################ RSiena ----
################################################################################ version September
################################################################################ 8, 2020 The
################################################################################ introductory
################################################################################ script is divided
################################################################################ into the following
################################################################################ script files:
################################################################################ Rscript01DataFormat.R,
################################################################################ followed by
################################################################################ RScriptSNADescriptives.R,
################################################################################ code for
################################################################################ descriptive
################################################################################ analysis of the
################################################################################ data, and
################################################################################ Rscript02SienaVariableFormat.R,
################################################################################ that formats data
################################################################################ and specifies the
################################################################################ model, and
################################################################################ Rscript03SienaRunModel.R,
################################################################################ that runs the
################################################################################ model and
################################################################################ estimates
################################################################################ parameters
################################################################################ Rscript04SienaBehaviour.R,
################################################################################ that illustrates
################################################################################ an example of
################################################################################ analysing the
################################################################################ coevolution of
################################################################################ networks and
################################################################################ behaviour Written
################################################################################ with contributions
################################################################################ by Robin Gauthier,
################################################################################ Tom Snijders, Ruth
################################################################################ Ripley, and Johan
################################################################################ Koskinen.
# Here is a short script for analysing the co-evolution of the friendship network and drinking
# behaviour for the 50 girls in the Teenage Friends and Lifestyle Study data
# (http://www.stats.ox.ac.uk/~snijders/siena/s50_data.zip), described in
# http://www.stats.ox.ac.uk/~snijders/siena/s50_data.htm
# Read in the adjacency matrices, covariates and dependent behavioural variable assuming data are
# in current working directory
# \t\tfriend.data.w1 <- as.matrix(read.table('s50-network1.dat')) # network \t\tfriend.data.w2 <-
# as.matrix(read.table('s50-network2.dat')) \t\tfriend.data.w3 <-
# as.matrix(read.table('s50-network3.dat')) \t\tdrink <- as.matrix(read.table('s50-alcohol.dat')) #
# behaviour \t\tsmoke <- as.matrix(read.table('s50-smoke.dat')) # covariate
# If you wish to make it easier, you can use this data set as included in the package - but the
# above is included to show you how to use data from files. To use the internal data set:
friend.data.w1 <- s501
friend.data.w2 <- s502
friend.data.w3 <- s503
drink <- s50a
smoke <- s50s
# At this point it is a good idea to use the sna package to plot the networks and the behavioural
# variable. Descriptive measures of the similarity of 'friends' with respect to behaviour (like
# Moran's I) are given by the function nacf() in the sna package. See the script
# RscriptSNADescriptives.R.
# Tell RSiena that the adjacency matrices are network data and in what order they should be treated
friendship <- sienaDependent(array(c(friend.data.w1, friend.data.w2, friend.data.w3), dim = c(50, 50,
3))) # create dependent variable
# Tell RSiena that the variable 'drink' should be treated as a dependent variable
drinkingbeh <- sienaDependent(drink, type = "behavior")
smoke1 <- coCovar(smoke[, 1])
# Define the data set and obtain the basic effects object
myCoEvolutionData <- sienaDataCreate(friendship, smoke1, drinkingbeh)
myCoEvolutionEff <- getEffects(myCoEvolutionData)
# Run reports to check that data is properly formated and to get some basic descriptives
print01Report(myCoEvolutionData, modelname = "s50_3_CoEvinit")
# Define the effects to include in the coevolution model Start with some structural effects (use
# the shortnames that you find in effectsDocumentation(myCoEvolutionEff) )
myCoEvolutionEff <- includeEffects(myCoEvolutionEff, transTrip, cycle3)#> effectName include fix test initialValue parm
#> 1 transitive triplets TRUE FALSE FALSE 0 0
#> 2 3-cycles TRUE FALSE FALSE 0 0# Include a homophily effect for the constant covariate smoking
myCoEvolutionEff <- includeEffects(myCoEvolutionEff, simX, interaction1 = "smoke1")#> effectName include fix test initialValue parm
#> 1 smoke1 similarity TRUE FALSE FALSE 0 0# If we want to parse out whether there is a selection or influence (or both) effect for drinking
# behaviour, we need to also include sender, receiver and homophily effects of drinking for
# friendship formation:
myCoEvolutionEff <- includeEffects(myCoEvolutionEff, egoX, altX, simX, interaction1 = "drinkingbeh")#> effectName include fix test initialValue parm
#> 1 drinkingbeh alter TRUE FALSE FALSE 0 0
#> 2 drinkingbeh ego TRUE FALSE FALSE 0 0
#> 3 drinkingbeh similarity TRUE FALSE FALSE 0 0# For the influence part, i.e. the effect of the network on behaviour, we specify the following
# effects: indegree, outdegree and assimilation effects for drinking
myCoEvolutionEff <- includeEffects(myCoEvolutionEff, name = "drinkingbeh", avAlt, indeg, outdeg, interaction1 = "friendship")#> effectName include fix test initialValue parm
#> 1 drinkingbeh indegree TRUE FALSE FALSE 0 0
#> 2 drinkingbeh outdegree TRUE FALSE FALSE 0 0
#> 3 drinkingbeh average alter TRUE FALSE FALSE 0 0# Check what effects you have decided to include:
myCoEvolutionEff#> name effectName include fix test initialValue parm
#> 1 friendship constant friendship rate (period 1) TRUE FALSE FALSE 4.69604 0
#> 2 friendship constant friendship rate (period 2) TRUE FALSE FALSE 4.32885 0
#> 3 friendship outdegree (density) TRUE FALSE FALSE -1.46770 0
#> 4 friendship reciprocity TRUE FALSE FALSE 0.00000 0
#> 5 friendship transitive triplets TRUE FALSE FALSE 0.00000 0
#> 6 friendship 3-cycles TRUE FALSE FALSE 0.00000 0
#> 7 friendship smoke1 similarity TRUE FALSE FALSE 0.00000 0
#> 8 friendship drinkingbeh alter TRUE FALSE FALSE 0.00000 0
#> 9 friendship drinkingbeh ego TRUE FALSE FALSE 0.00000 0
#> 10 friendship drinkingbeh similarity TRUE FALSE FALSE 0.00000 0
#> 11 drinkingbeh rate drinkingbeh (period 1) TRUE FALSE FALSE 0.70571 0
#> 12 drinkingbeh rate drinkingbeh (period 2) TRUE FALSE FALSE 0.84939 0
#> 13 drinkingbeh drinkingbeh linear shape TRUE FALSE FALSE 0.32237 0
#> 14 drinkingbeh drinkingbeh quadratic shape TRUE FALSE FALSE 0.00000 0
#> 15 drinkingbeh drinkingbeh indegree TRUE FALSE FALSE 0.00000 0
#> 16 drinkingbeh drinkingbeh outdegree TRUE FALSE FALSE 0.00000 0
#> 17 drinkingbeh drinkingbeh average alter TRUE FALSE FALSE 0.00000 0# Now we have to define the algorithm settings. The defaults are adequate. You only have to
# specify the filename that will receive the results in text format.
myCoEvAlgorithm <- sienaAlgorithmCreate(projname = "s50CoEv_3")#> If you use this algorithm object, siena07 will create/use an output file s50CoEv_3.txt .# Finally, estimate the model; the whole command is put in parentheses to have the results printed
# directly to the screen.
(ans <- siena07(myCoEvAlgorithm, data = myCoEvolutionData, effects = myCoEvolutionEff))#> Estimates, standard errors and convergence t-ratios
#>
#> Estimate Standard Convergence
#> Error t-ratio
#> Network Dynamics
#> 1. rate constant friendship rate (period 1) 6.4800 ( 1.2439 ) -0.0641
#> 2. rate constant friendship rate (period 2) 5.1314 ( 0.8404 ) 0.0138
#> 3. eval outdegree (density) -2.7681 ( 0.1556 ) 0.0127
#> 4. eval reciprocity 2.3747 ( 0.2222 ) -0.0011
#> 5. eval transitive triplets 0.6483 ( 0.1344 ) 0.0402
#> 6. eval 3-cycles -0.0730 ( 0.2837 ) 0.0360
#> 7. eval smoke1 similarity 0.1858 ( 0.2230 ) 0.0932
#> 8. eval drinkingbeh alter -0.0383 ( 0.1164 ) -0.0521
#> 9. eval drinkingbeh ego 0.0720 ( 0.1265 ) -0.0709
#> 10. eval drinkingbeh similarity 1.3180 ( 0.6331 ) 0.0602
#>
#> Behavior Dynamics
#> 11. rate rate drinkingbeh (period 1) 1.2443 ( 0.3677 ) 0.0144
#> 12. rate rate drinkingbeh (period 2) 1.6758 ( 0.4256 ) -0.1003
#> 13. eval drinkingbeh linear shape -0.8497 ( 2.2964 ) 0.0008
#> 14. eval drinkingbeh quadratic shape -1.7646 ( 3.7417 ) -0.0192
#> 15. eval drinkingbeh indegree -0.9942 ( 2.8822 ) 0.0246
#> 16. eval drinkingbeh outdegree 1.7308 ( 4.3426 ) 0.0272
#> 17. eval drinkingbeh average alter 3.8412 ( 8.6979 ) 0.0023
#>
#> Overall maximum convergence ratio: 0.2272
#>
#>
#> Total of 3406 iteration steps.# THE RESULTS
# Note that the 'convergence t-ratio' is the t-ratio for convergence checking, not the t statistic
# for testing the significance of this effect. (See Section 6.1.2 of the manual.) For good
# convergence, the t-ratios for convergence all should be less than .1 in absolute value, and the
# overall maximum convergence ratio should be less than 0.25. If this is not yet the case, you
# should try again, starting from the last estimate as the 'previous answer':
(ans1 <- siena07(myCoEvAlgorithm, data = myCoEvolutionData, effects = myCoEvolutionEff, prevAns = ans))#> Estimates, standard errors and convergence t-ratios
#>
#> Estimate Standard Convergence
#> Error t-ratio
#> Network Dynamics
#> 1. rate constant friendship rate (period 1) 6.5731 ( 1.3339 ) 0.0203
#> 2. rate constant friendship rate (period 2) 5.1349 ( 1.0727 ) -0.0655
#> 3. eval outdegree (density) -2.7726 ( 0.1761 ) -0.0059
#> 4. eval reciprocity 2.3789 ( 0.2463 ) -0.0020
#> 5. eval transitive triplets 0.6245 ( 0.1544 ) -0.0103
#> 6. eval 3-cycles -0.0286 ( 0.3158 ) 0.0065
#> 7. eval smoke1 similarity 0.1797 ( 0.2613 ) 0.0466
#> 8. eval drinkingbeh alter -0.0378 ( 0.1149 ) 0.0249
#> 9. eval drinkingbeh ego 0.0744 ( 0.1153 ) 0.0125
#> 10. eval drinkingbeh similarity 1.3572 ( 0.8138 ) 0.0581
#>
#> Behavior Dynamics
#> 11. rate rate drinkingbeh (period 1) 1.2365 ( 0.3339 ) 0.0308
#> 12. rate rate drinkingbeh (period 2) 1.7081 ( 0.4566 ) 0.0222
#> 13. eval drinkingbeh linear shape -0.9161 ( 2.0818 ) -0.0003
#> 14. eval drinkingbeh quadratic shape -1.7939 ( 2.8311 ) -0.0061
#> 15. eval drinkingbeh indegree -1.0705 ( 2.4082 ) 0.0511
#> 16. eval drinkingbeh outdegree 1.8347 ( 3.6174 ) 0.0336
#> 17. eval drinkingbeh average alter 3.8696 ( 6.2498 ) -0.0410
#>
#> Overall maximum convergence ratio: 0.2656
#>
#>
#> Total of 3348 iteration steps.# which can be repeated if necessary.
# For this small data set, the model for behavior dynamics is over-specified, leading to some very
# large standard errors. For this data set it is better to drop the degree effects on behaviour,
# because the data does not contain enough information to estimate them.
myCoEvolutionEff2 <- includeEffects(myCoEvolutionEff, name = "drinkingbeh", indeg, outdeg, interaction1 = "friendship",
include = FALSE)#> [1] effectName include fix test initialValue parm
#> <0 rows> (or 0-length row.names)(ans2 <- siena07(myCoEvAlgorithm, data = myCoEvolutionData, effects = myCoEvolutionEff2))#> Estimates, standard errors and convergence t-ratios
#>
#> Estimate Standard Convergence
#> Error t-ratio
#> Network Dynamics
#> 1. rate constant friendship rate (period 1) 6.5490 ( 1.1566 ) -0.0211
#> 2. rate constant friendship rate (period 2) 5.1581 ( 0.9815 ) 0.0637
#> 3. eval outdegree (density) -2.7601 ( 0.1450 ) 0.0188
#> 4. eval reciprocity 2.3730 ( 0.2064 ) 0.0096
#> 5. eval transitive triplets 0.6382 ( 0.1406 ) 0.0062
#> 6. eval 3-cycles -0.0604 ( 0.2859 ) 0.0116
#> 7. eval smoke1 similarity 0.1764 ( 0.2209 ) 0.0110
#> 8. eval drinkingbeh alter -0.0436 ( 0.1213 ) -0.0188
#> 9. eval drinkingbeh ego 0.0717 ( 0.1337 ) -0.0115
#> 10. eval drinkingbeh similarity 1.3418 ( 0.6484 ) -0.0177
#>
#> Behavior Dynamics
#> 11. rate rate drinkingbeh (period 1) 1.3194 ( 0.3736 ) -0.0149
#> 12. rate rate drinkingbeh (period 2) 1.7806 ( 0.5254 ) -0.0173
#> 13. eval drinkingbeh linear shape 0.4179 ( 0.2382 ) -0.0202
#> 14. eval drinkingbeh quadratic shape -0.5948 ( 0.3620 ) 0.0038
#> 15. eval drinkingbeh average alter 1.3014 ( 0.9056 ) -0.0222
#>
#> Overall maximum convergence ratio: 0.1123
#>
#>
#> Total of 3252 iteration steps.############################################################################### Some other effects
############################################################################### ## The set of
############################################################################### available effects
############################################################################### for this data set
############################################################################### can be obtained by
############################################################################### requesting
effectsDocumentation(myCoEvolutionEff)
# See the manual, Chapter 12, for the meaning of these effects. To study the direct effect of the
# actor covariate smoking on the dependent variable drinking, use the effFrom effect:
myCoEvolutionEff3 <- includeEffects(myCoEvolutionEff2, name = "drinkingbeh", effFrom, interaction1 = "smoke1")#> effectName include fix test initialValue parm
#> 1 drinkingbeh: effect from smoke1 TRUE FALSE FALSE 0 0# Since we already have a good result for a simpler model, it is helpful to start estimating from
# these estimates as starting values:
(ans3 <- siena07(myCoEvAlgorithm, data = myCoEvolutionData, effects = myCoEvolutionEff3, prevAns = ans2))#> Estimates, standard errors and convergence t-ratios
#>
#> Estimate Standard Convergence
#> Error t-ratio
#> Network Dynamics
#> 1. rate constant friendship rate (period 1) 6.5686 ( 1.2079 ) -0.0390
#> 2. rate constant friendship rate (period 2) 5.1513 ( 0.8772 ) -0.0364
#> 3. eval outdegree (density) -2.7731 ( 0.1439 ) -0.0535
#> 4. eval reciprocity 2.3938 ( 0.2113 ) -0.0156
#> 5. eval transitive triplets 0.6305 ( 0.1458 ) -0.0261
#> 6. eval 3-cycles -0.0424 ( 0.2785 ) -0.0062
#> 7. eval smoke1 similarity 0.1868 ( 0.2309 ) 0.0027
#> 8. eval drinkingbeh alter -0.0394 ( 0.1231 ) -0.0175
#> 9. eval drinkingbeh ego 0.0679 ( 0.1239 ) -0.0448
#> 10. eval drinkingbeh similarity 1.3310 ( 0.6729 ) 0.0033
#>
#> Behavior Dynamics
#> 11. rate rate drinkingbeh (period 1) 1.2809 ( 0.3422 ) -0.0173
#> 12. rate rate drinkingbeh (period 2) 1.7066 ( 0.4539 ) 0.0096
#> 13. eval drinkingbeh linear shape 0.4345 ( 0.2945 ) -0.0118
#> 14. eval drinkingbeh quadratic shape -0.6347 ( 0.4529 ) 0.0359
#> 15. eval drinkingbeh average alter 1.5132 ( 1.1349 ) 0.0848
#> 16. eval drinkingbeh: effect from smoke1 -0.2412 ( 0.3955 ) -0.0025
#>
#> Overall maximum convergence ratio: 0.2508
#>
#>
#> Total of 3322 iteration steps.# In my case, convergence was not good enough:
(ans3 <- siena07(myCoEvAlgorithm, data = myCoEvolutionData, effects = myCoEvolutionEff3, prevAns = ans3))#> Estimates, standard errors and convergence t-ratios
#>
#> Estimate Standard Convergence
#> Error t-ratio
#> Network Dynamics
#> 1. rate constant friendship rate (period 1) 6.5114 ( 1.0970 ) -0.0251
#> 2. rate constant friendship rate (period 2) 5.2016 ( 0.8726 ) 0.0756
#> 3. eval outdegree (density) -2.7703 ( 0.1311 ) -0.0072
#> 4. eval reciprocity 2.3784 ( 0.2060 ) 0.0108
#> 5. eval transitive triplets 0.6279 ( 0.1449 ) 0.0256
#> 6. eval 3-cycles -0.0341 ( 0.2875 ) 0.0401
#> 7. eval smoke1 similarity 0.1822 ( 0.2057 ) -0.0191
#> 8. eval drinkingbeh alter -0.0414 ( 0.1202 ) 0.0656
#> 9. eval drinkingbeh ego 0.0744 ( 0.1268 ) 0.0304
#> 10. eval drinkingbeh similarity 1.3608 ( 0.6272 ) 0.0243
#>
#> Behavior Dynamics
#> 11. rate rate drinkingbeh (period 1) 1.2897 ( 0.3725 ) 0.0816
#> 12. rate rate drinkingbeh (period 2) 1.7271 ( 0.4554 ) 0.0116
#> 13. eval drinkingbeh linear shape 0.4472 ( 0.2380 ) 0.0261
#> 14. eval drinkingbeh quadratic shape -0.6494 ( 0.3851 ) -0.0409
#> 15. eval drinkingbeh average alter 1.5168 ( 0.9566 ) 0.0366
#> 16. eval drinkingbeh: effect from smoke1 -0.2083 ( 0.3921 ) 0.0309
#>
#> Overall maximum convergence ratio: 0.2431
#>
#>
#> Total of 3267 iteration steps.# You can get a nicer presentation of the results in a file in your working directory in LaTeX by
siena.table(ans3)
# and in html (can be imported into MS-Word) by
siena.table(ans3, type = "html")Copyright © 2020 Jochem Tolsma