# Chapter 6amen packages

## 6.1 Outline

• Gaussian AME model: ame
• Different relation: ame(...,model=,...)
• ordinal data
• censored and fixed rank nomination data
• sampled or missing data
• symmetric relation: ame(...,symmetric=TRUE,...)
• repeated measures data: longitudal data ame_rep(Y,Xdyad,Xrow,Xcol)

## 6.2ame()

$z_{i,j}=\beta_d^Tx_{d,i,j}+\beta_r^Tx_{r,i}+\beta_c^Tx_{c,j}+\mu+a_i+b_j+u_i^Tv_j+\epsilon_{i,j}$

$y_{i,j}=g(z_{i,j})$

ame(Y,Xd,Xr,Xc,model,R,rvar=TRUE,cvar=TRUE,dcor=TRUE,symmetric=FALSE)

model:

• “nrm”: continuous Y
• “bin”: binary Y
• “ord”: ordinal Y
• “cbin”: censored binary Y
• “frn”: fixed rank Y
• “rrl”: row ranks (relationships across rows of Y are not directly comparable)

• Y: named matrix $$n\times n \times p_d$$
• Xd: named array (dyadic covariates)
• Xr: named matrix (nodal covariates)
• Xc: named matrix (nodal covariates)
• R: number (dimension of latent factor)

• rvar,cvar,dcor: whether to include row/column/corr effect or not
• symmetric: symmetric outcome or not

## 6.3 Input: prepare the dataset

prepare the dataset

• Get the Xnode,Xdyad from igraph object.
• Get the Xnode,Xdyad from edgelist dataframe.
library(igraph)
library(igraphdata)
library(dplyr)

From igraph object:

data("USairports")
?USairports
USairports
## IGRAPH bf6202d DN-- 755 23473 -- US airports
## + attr: name (g/c), name (v/c), City (v/c), Position (v/c),
## | Carrier (e/c), Departures (e/n), Seats (e/n), Passengers (e/n),
## | Aircraft (e/n), Distance (e/n)
## + edges from bf6202d (vertex names):
##  [1] BGR->JFK BGR->JFK BOS->EWR ANC->JFK JFK->ANC LAS->LAX MIA->JFK
##  [8] EWR->ANC BJC->MIA MIA->BJC TEB->ANC JFK->LAX LAX->JFK LAX->SFO
## [15] AEX->LAS BFI->SBA ELM->PIT GEG->SUN ICT->PBI LAS->LAX LAS->PBI
## [22] LAS->SFO LAX->LAS PBI->AEX PBI->ICT PIT->VCT SFO->LAX VCT->DWH
## [29] IAD->JFK ABE->CLT ABE->HPN AGS->CLT AGS->CLT AVL->CLT AVL->CLT
## [36] AVP->CLT AVP->PHL BDL->CLT BHM->CLT BHM->CLT BNA->CLT BNA->CLT
## + ... omitted several edges
#For vertex attributes -Xnode
Xnode_ls=vertex_attr(USairports)
Xnode=matrix(unlist(Xnode_ls),ncol=length(Xnode_ls))
colnames(Xnode)=names(Xnode_ls)
rownames(Xnode)=Xnode_ls$name Xnode=Xnode[,-1] #For edge attributes - Xdyad Xdyad_ls=edge_attr(USairports) Xdyad=array(c(as_adjacency_matrix(USairports,sparse = FALSE),unlist(Xdyad_ls)),dim=c(nrow(Xnode),nrow(Xnode),length(Xdyad_ls)+1)) dimnames(Xdyad)[[1]]=dimnames(Xdyad)[[2]]=Xnode_ls$name
dimnames(Xdyad)[[3]]=c("relation",names(Xdyad_ls))

From dataframe:

transfer it to a igraph then to the array may be the easiest.

#Xdyad
df=igraph::as_data_frame(USairports)
net=graph_from_data_frame(df)

vnames=V(net)$name Xdyad_ls=edge_attr(USairports) Xdyad=array(c(as_adjacency_matrix(USairports,sparse = FALSE),unlist(Xdyad_ls)),dim=c(nrow(Xnode),nrow(Xnode),length(Xdyad_ls)+1)) dimnames(Xdyad)[[1]]=dimnames(Xdyad)[[2]]=vnames dimnames(Xdyad)[[3]]=c("relation",names(Xdyad_ls)) #Xnode #If your nodal attributes stored in a data.frame Xnode_df=as.data.frame(Xnode)%>%mutate(Vertex=rownames(Xnode))%>%arrange(Vertex) head(Xnode_df) ## City Position Vertex ## 1 Peach Springs, AZ N355925 W1134859 1G4 ## 2 Bradley Lake, AK N565311 W1340930 A23 ## 3 Pogo Mines, AK N592603 W1514228 A27 ## 4 Kiluda Bay, AK N570308 W1352046 A29 ## 5 Allentown/Bethlehem/Easton, PA N403909 W0752625 ABE ## 6 Abilene, TX N322441 W0994055 ABI #Xnode #need to match the order of Xnode to the order of the Xdyad Xnode_df=Xnode_df%>%arrange(match(Vertex,vnames)) Xnode=as.matrix(Xnode_df) rownames(Xnode)=Xnode[,"Vertex"] Xnode=Xnode[,-3] Xnode%>%head() ## City Position ## BGR "Bangor, ME" "N444827 W0684941" ## BOS "Boston, MA" "N422152 W0710019" ## ANC "Anchorage, AK" "N611028 W1495947" ## JFK "New York, NY" "N403823 W0734644" ## LAS "Las Vegas, NV" "N360449 W1150908" ## MIA "Miami, FL" "N254736 W0801726" ## 6.4 Fit the model: the Gaussian AME model (continuous Y) ### 6.4.1 The Gaussian AME model #install.packages("amen") library(amen) library(ggplot2) $$y_{i,j}=\mu+a_i+b_j+\epsilon_{i,j}$$ # use the trade data (export) from top 30 countries ranked by gdp data(IR90s) names(IR90s) # list ## [1] "dyadvars" "nodevars" dim(IR90s$dyadvars) #- array
## [1] 130 130   5
dimnames(IR90s$dyadvars) ## [[1]] ## [1] "AFG" "ALB" "ALG" "ANG" "ARG" "AUL" "AUS" "BAH" "BEL" "BEN" "BFO" ## [12] "BHU" "BNG" "BOL" "BOT" "BRA" "BUI" "BUL" "CAM" "CAN" "CAO" "CDI" ## [23] "CEN" "CHA" "CHL" "CHN" "COL" "COM" "CON" "COS" "CUB" "CYP" "DEN" ## [34] "DJI" "DOM" "DRC" "ECU" "EGY" "EQG" "FIN" "FJI" "FRN" "GAB" "GAM" ## [45] "GHA" "GNB" "GRC" "GUA" "GUI" "GUY" "HAI" "HON" "HUN" "IND" "INS" ## [56] "IRE" "IRN" "IRQ" "ISR" "ITA" "JAM" "JOR" "JPN" "KEN" "LAO" "LBR" ## [67] "LES" "LIB" "MAA" "MAG" "MAL" "MAS" "MAW" "MEX" "MLI" "MON" "MOR" ## [78] "MYA" "MZM" "NAM" "NEP" "NEW" "NIC" "NIG" "NIR" "NOR" "NTH" "OMA" ## [89] "PAK" "PAN" "PAR" "PHI" "PNG" "POL" "POR" "PRK" "QAT" "ROK" "RUM" ## [100] "RWA" "SAF" "SAL" "SAU" "SEN" "SIE" "SIN" "SOM" "SPN" "SRI" "SUD" ## [111] "SWA" "SWD" "SWZ" "SYR" "TAW" "TAZ" "THI" "TOG" "TRI" "TUN" "TUR" ## [122] "UAE" "UGA" "UKG" "URU" "USA" "VEN" "YEM" "ZAM" "ZIM" ## ## [[2]] ## [1] "AFG" "ALB" "ALG" "ANG" "ARG" "AUL" "AUS" "BAH" "BEL" "BEN" "BFO" ## [12] "BHU" "BNG" "BOL" "BOT" "BRA" "BUI" "BUL" "CAM" "CAN" "CAO" "CDI" ## [23] "CEN" "CHA" "CHL" "CHN" "COL" "COM" "CON" "COS" "CUB" "CYP" "DEN" ## [34] "DJI" "DOM" "DRC" "ECU" "EGY" "EQG" "FIN" "FJI" "FRN" "GAB" "GAM" ## [45] "GHA" "GNB" "GRC" "GUA" "GUI" "GUY" "HAI" "HON" "HUN" "IND" "INS" ## [56] "IRE" "IRN" "IRQ" "ISR" "ITA" "JAM" "JOR" "JPN" "KEN" "LAO" "LBR" ## [67] "LES" "LIB" "MAA" "MAG" "MAL" "MAS" "MAW" "MEX" "MLI" "MON" "MOR" ## [78] "MYA" "MZM" "NAM" "NEP" "NEW" "NIC" "NIG" "NIR" "NOR" "NTH" "OMA" ## [89] "PAK" "PAN" "PAR" "PHI" "PNG" "POL" "POR" "PRK" "QAT" "ROK" "RUM" ## [100] "RWA" "SAF" "SAL" "SAU" "SEN" "SIE" "SIN" "SOM" "SPN" "SRI" "SUD" ## [111] "SWA" "SWD" "SWZ" "SYR" "TAW" "TAZ" "THI" "TOG" "TRI" "TUN" "TUR" ## [122] "UAE" "UGA" "UKG" "URU" "USA" "VEN" "YEM" "ZAM" "ZIM" ## ## [[3]] ## [1] "conflicts" "exports" "distance" "shared_igos" "polity_int" gdp=IR90s$nodevars[,2]
topgdp=which(gdp>=sort(gdp,decreasing = TRUE)[30])
Y=log(IR90s$dyadvars[topgdp,topgdp,2]+1) Y[1:5,1:5] ## ARG AUL BEL BNG BRA ## ARG NA 0.05826891 0.2468601 0.03922071 1.76473080 ## AUL 0.0861777 NA 0.3784364 0.10436002 0.21511138 ## BEL 0.2700271 0.35065687 NA 0.01980263 0.39877612 ## BNG 0.0000000 0.01980263 0.1222176 NA 0.01980263 ## BRA 1.6937791 0.23901690 0.6205765 0.03922071 NA ### 6.4.2 Social relations model (SRM) $y_{i,j}=\mu+a_i+b_j+\epsilon_{i,j}$ Input: Y - a named matrix # fit the model fit_SRM=ame(Y,model="nrm",plot=FALSE,print = FALSE) #normal AME model; by default: niter=10,000; save every 25; burnin=500 #output names(fit_SRM) ## [1] "BETA" "VC" "APM" "BPM" "U" "V" "UVPM" "EZ" "YPM" "GOF" mean(fit_SRM$BETA) #mu_hat
## [1] 0.6616995
(muhat=mean(Y,na.rm = TRUE))
## [1] 0.680044
apply(fit_SRM$VC,2,mean) # covariance  ## va cab vb rho ve ## 0.2809375 0.2364107 0.2722367 0.8577201 0.2387065 ahat=rowMeans(Y,na.rm = TRUE)-mean(Y,na.rm = TRUE) bhat=colMeans(Y,na.rm = TRUE)-mean(Y,na.rm = TRUE) (corab=cov(cbind(ahat,bhat))) ## ahat bhat ## ahat 0.2407563 0.2290788 ## bhat 0.2290788 0.2289489 (R=Y-(muhat+outer(ahat,bhat,"+"))) ## ARG AUL BEL BNG BRA ## ARG NA -0.138060630 -0.23544663 0.312875331 1.63527007 ## AUL -0.16956545 NA -0.43273546 0.049149445 -0.24321454 ## BEL -0.21364558 -0.402467428 NA -0.263337516 -0.28747937 ## BNG 0.28580442 0.036155468 -0.14740670 NA 0.10302428 ## BRA 1.46096778 -0.263245960 -0.16766355 0.006942009 NA ## CAN -0.17367495 0.008106998 -0.30453842 -0.086801691 0.10095129 ## CHN -0.28802510 0.388372858 -0.30185657 -0.100071302 -0.36177499 ## COL 0.29564048 -0.022129678 -0.17829524 0.437904147 0.14811749 ## EGY 0.25762273 -0.001878520 -0.21162833 0.468105637 0.07484259 ## FRN -0.23537558 -0.357291724 1.60148082 -0.635013732 -0.29042291 ## IND -0.01417791 -0.043803438 0.06583220 0.637489812 -0.13599932 ## INS -0.13556799 0.325168886 -0.32399791 0.102005862 -0.26539140 ## IRN 0.03803907 -0.221462181 -0.21728509 0.248521975 0.20701985 ## ITA -0.04846778 -0.160389153 0.55599067 -0.545832983 0.13043260 ## JPN -0.77132395 0.587375094 -0.43066468 -0.739867701 -0.32622191 ## MEX 0.31445810 -0.176365556 -0.33770028 0.171400968 0.25071177 ## NTH -0.35738098 -0.398742583 2.11143771 -0.261291065 -0.33187377 ## PAK 0.23763171 0.043687688 -0.21538203 0.504539360 0.04480512 ## PHI 0.06168391 -0.043764695 -0.40838500 0.252266154 -0.11162848 ## POL 0.20701146 -0.072196296 0.01406599 0.387935563 0.10589815 ## ROK -0.16700865 0.192608090 -0.72935785 -0.075675369 -0.15336832 ## SAF 0.15365832 -0.024443477 0.02378144 0.214428959 0.01016447 ## SAU -0.31498022 -0.149696308 -0.66524792 -0.029087793 0.26547920 ## SPN 0.27723450 -0.300832361 0.37460163 -0.073409842 0.03430044 ## SWD -0.10048691 0.088730653 0.60848888 -0.069246392 -0.02568636 ## TAW -0.26831811 0.249196248 -0.55259409 -0.035524670 -0.25764209 ## THI -0.10070850 0.175886982 -0.21970992 0.118516213 -0.24823149 ## TUR 0.13431957 -0.106144478 -0.09025494 0.324999845 -0.04884512 ## UKG -0.58169724 0.242724120 1.08346413 -0.559239933 -0.46004035 ## USA -0.16801641 0.489417309 -0.09654737 -1.359712929 0.35838719 ## CAN CHN COL EGY FRN ## ARG -0.25750202 -0.218112259 0.379283330 0.335425500 -0.50910133 ## AUL -0.14391354 0.487239680 -0.095293669 0.110238741 -0.63291671 ## BEL -0.54194908 -0.388986507 -0.266064828 -0.196471710 1.70626150 ## BNG -0.11692294 -0.394894646 0.426451271 0.384104903 -0.44483896 ## BRA -0.12658392 -0.164793316 0.229919068 0.037860427 -0.12880792 ## CAN NA 0.192689487 0.053983278 -0.198984116 -0.33600377 ## CHN 0.38302118 NA -0.495164243 -0.354867835 0.13190355 ## COL -0.06750000 -0.482167329 NA 0.325721722 -0.50130398 ## EGY -0.19296065 -0.432645788 0.398269579 NA -0.45714730 ## FRN -0.10102537 -0.153364961 -0.500193843 -0.216018305 NA ## IND -0.19454812 0.059967577 0.088369094 0.172028032 -0.36779412 ## INS -0.30654846 0.403748876 -0.033387417 0.016730889 -0.53243819 ## IRN -0.31901826 -0.319886607 0.178685918 0.136339550 -0.15494203 ## ITA -0.16692944 0.001411349 -0.407648796 -0.070228744 1.55611087 ## JPN 0.46908045 1.458427987 -0.647484216 -0.878718373 -0.02644183 ## MEX 1.10390222 -0.527509119 0.440818988 0.069168873 -0.53184292 ## NTH -0.50556372 -0.410390081 -0.312044308 -0.159546595 1.24932142 ## PAK -0.12152866 -0.112884023 0.358475938 0.325981866 -0.45832611 ## PHI -0.05400932 0.025341858 0.192380427 0.150034059 -0.62442162 ## POL -0.20675566 -0.438707890 0.328049837 0.333461456 -0.07819955 ## ROK 0.12887478 1.397776911 -0.249505355 -0.227078508 -0.54070127 ## SAF -0.24137825 -0.190864320 0.184102035 0.112196864 -0.42706692 ## SAU -0.43764211 -0.523791529 -0.184283699 0.056089216 0.02690898 ## SPN -0.36535493 -0.457920467 0.062013184 -0.006544036 1.60143892 ## SWD -0.14378535 -0.124608851 -0.045362235 0.028029281 0.30922160 ## TAW 0.17100764 1.413411595 -0.240344421 -0.213304377 -0.48486254 ## THI -0.09287004 0.443899383 -0.008380224 -0.012358346 -0.32430698 ## TUR -0.25150906 -0.416121434 0.265114119 0.417978469 -0.01949016 ## UKG 0.30907514 -0.173171638 -0.577625928 -0.451674396 1.21014874 ## USA 2.44390557 0.493707570 -0.071731672 -0.594444456 0.38161638 ## IND INS IRN ITA JPN ## ARG 0.12233225 0.062477771 0.503935683 -0.207345734 -0.720363646 ## AUL 0.31608726 0.612026831 0.111444799 -0.377986917 1.199289037 ## BEL 0.54742550 -0.369674544 -0.087705809 0.972810738 -0.706327549 ## BNG 0.20872090 0.162892849 0.377402340 -0.347264973 -0.831081422 ## BRA -0.11009861 -0.126656283 0.148869672 0.005514737 -0.005172175 ## CAN -0.21823527 -0.065374795 0.010165397 -0.260472972 0.698722565 ## CHN -0.16138315 -0.053643104 -0.334087963 -0.008877084 1.527796818 ## COL 0.13101767 0.104509668 0.289748777 -0.412091633 -0.697192715 ## EGY 0.31527181 0.134711157 0.309999935 0.149907060 -0.850213278 ## FRN -0.45419538 -0.482181747 -0.291211963 1.563964415 -0.159540500 ## IND NA 0.121773776 0.136055086 -0.179982776 -0.061472319 ## INS 0.09532398 NA -0.064498647 -0.446133793 1.259109783 ## IRN 0.39796902 0.136288886 NA 0.157305726 0.248951379 ## ITA -0.33484659 -0.451783448 0.022689324 NA -0.186033475 ## JPN -0.37347476 0.522930506 -0.460865139 -0.490816315 NA ## MEX -0.09664567 -0.103724604 0.043048224 -0.670364313 -0.285816186 ## NTH -0.30663112 -0.268753590 -0.188139432 1.146405249 -0.763248726 ## PAK 0.16918351 0.189327200 0.328475202 -0.407143096 -0.602324847 ## PHI -0.03501185 -0.004167284 0.226328416 -0.631346237 0.408153473 ## POL 0.11941733 0.074149555 0.269050575 0.134876289 -0.930383351 ## ROK -0.04511087 0.403568310 -0.003768855 -0.581148577 0.979043460 ## SAF 0.08821251 -0.023655341 0.115063752 0.397588251 0.001665143 ## SAU 0.60958685 0.072835799 -0.141525081 -0.135394695 0.862585232 ## SPN -0.26162927 -0.235986589 0.013274379 1.037508198 -0.735832846 ## SWD -0.20862783 -0.154567937 -0.031849185 0.184997624 -0.406276962 ## TAW -0.23988349 0.254384458 -0.328614065 -0.532525363 0.897056769 ## THI -0.09296008 0.208935980 0.050403550 -0.480267175 0.822712249 ## TUR 0.06587136 0.020875748 0.459563758 0.356238545 -0.915357528 ## UKG 0.03681840 -0.492835042 -0.348181213 0.784074660 -0.190145269 ## USA -0.51389508 -0.366074511 -1.484337802 -0.238101334 1.189580283 ## MEX NTH PAK PHI POL ## ARG 0.124973579 -0.163851688 0.26677001 0.042982916 0.119797247 ## AUL -0.310989745 -0.615849561 0.10277947 0.207138727 -0.228299303 ## BEL -0.401049524 1.634785609 -0.23496497 -0.393852530 0.068097435 ## BNG 0.090795880 -0.306103608 0.46040420 0.245812999 0.293639793 ## BRA 0.275277690 -0.001306934 -0.02977357 -0.092132674 -0.077922490 ## CAN 0.452297063 -0.454111393 -0.14273855 -0.223616702 -0.253081581 ## CHN -0.269940422 0.057863812 -0.23060758 -0.080421871 -0.271563669 ## COL 0.164580663 -0.295447541 0.37275063 0.187429818 0.264815414 ## EGY 0.062614189 -0.360718556 0.41280442 0.207680976 0.255507771 ## FRN -0.318888621 0.799615236 -0.53671417 -0.645662010 0.084221278 ## IND -0.145816836 -0.343502900 0.15790252 0.050543246 0.047076746 ## INS -0.181831266 -0.237372227 0.06375980 0.154997473 -0.082623168 ## IRN -0.137166846 -0.051777017 0.33893257 0.288201907 0.192927858 ## ITA -0.248226920 0.258244786 -0.38974953 -0.671144170 0.491961966 ## JPN 0.004025465 -0.109779114 -0.52063404 0.410646473 -0.916539749 ## MEX NA -0.513941502 0.11619778 -0.059270735 -0.011443940 ## NTH -0.511898158 NA -0.26420227 -0.351083013 0.238341064 ## PAK 0.042623174 -0.365420878 NA 0.197446137 0.225664460 ## PHI 0.022239475 -0.038082507 0.18711264 NA 0.069421246 ## POL 0.031234280 0.071721116 0.33263434 0.157069705 NA ## ROK 0.129309204 -0.547536198 -0.11801421 0.327776769 -0.121976834 ## SAF -0.132321994 -0.353130475 0.17883425 0.003175342 0.060571588 ## SAU -0.471148926 0.055256995 0.30514475 0.313262336 -0.307242881 ## SPN 0.194343258 0.360879345 -0.10009708 -0.276028150 0.186566145 ## SWD -0.179535571 0.503706625 -0.08699767 -0.253969344 0.308312843 ## TAW -0.052411936 0.187455423 -0.18336990 0.358358017 -0.270433067 ## THI -0.138874566 -0.093140180 0.08973034 0.169517278 -0.074914667 ## TUR -0.050932801 -0.069031625 0.31755465 0.084377812 0.216761995 ## UKG -0.554981679 1.055749094 -0.41147145 -0.638756464 0.008691453 ## USA 2.229606283 0.171599983 -0.98303369 -0.083441062 -1.005014825 ## ROK SAF SAU SPN SWD ## ARG -0.362233457 0.19553986 3.860386e-05 0.14735374 -0.06359310 ## AUL 0.860561467 0.07592911 -5.323656e-02 -0.54508562 -0.30065074 ## BEL -0.598980861 -0.14503179 -2.107145e-01 0.58202663 0.44895404 ## BNG -0.270629009 0.25980223 1.834506e-01 -0.21711215 0.14908928 ## BRA -0.121161499 -0.05996261 1.255901e-02 -0.05215012 -0.20173297 ## CAN 0.181178447 -0.19144664 -1.989677e-01 -0.50389239 -0.30121427 ## CHN 0.727143852 -0.32029690 -2.660169e-01 -0.11278727 -0.22101478 ## COL -0.289791477 0.21136938 1.448700e-01 -0.14038186 0.13839675 ## EGY -0.186543399 0.26117934 3.441694e-01 -0.10346365 0.15864791 ## FRN -0.476336945 -0.28735326 -2.090332e-01 1.54531048 0.35832183 ## IND -0.135396692 0.06611456 2.441142e-01 -0.25510686 -0.10935241 ## INS 0.674863870 -0.13366139 8.669695e-03 -0.20321934 -0.29281171 ## IRN 0.313222493 0.46946173 -2.547492e-02 -0.03124152 0.07721459 ## ITA -0.349938573 -0.19521948 -7.890806e-02 1.08235528 0.12418180 ## JPN 1.312006916 -0.26018320 -8.512040e-02 -0.36112944 -0.39401347 ## MEX -0.356498469 -0.01591308 -9.264560e-02 0.17996748 -0.17550867 ## NTH -0.477808270 -0.14222111 -2.415952e-01 0.62127121 0.81077717 ## PAK -0.141540504 0.21162953 2.452869e-01 -0.19411571 0.10955188 ## PHI -0.007832581 0.07452155 -1.178041e-02 -0.40466298 -0.08498157 ## POL -0.321272456 0.16140080 1.043692e-01 -0.21725033 0.44620262 ## ROK NA -0.26446981 -3.569072e-02 -0.39713177 -0.40702337 ## SAF 0.115341994 NA 4.135417e-02 -0.12570070 -0.12281876 ## SAU 0.907089702 -0.02884924 NA -0.13469994 -0.34386266 ## SPN -0.530825059 -0.12089638 -7.674545e-02 NA 0.12915560 ## SWD -0.421035710 -0.07847988 -8.955842e-02 0.11593934 NA ## TAW 0.113849725 -0.01693698 -2.010367e-01 -0.47075430 -0.27846171 ## THI -0.088784337 -0.03655959 -2.442148e-02 -0.28541641 -0.16095996 ## TUR -0.312104881 0.13730491 3.096975e-01 -0.05875336 0.06960934 ## UKG -0.383108089 0.28418233 2.598720e-01 0.77849980 0.75819531 ## USA 0.853908792 -0.43607760 1.619809e-02 -0.25615041 -0.33218743 ## TAW THI TUR UKG USA ## ARG -0.26786838 -0.016537857 0.07753578 -0.798550230 -0.73684143 ## AUL 0.53918578 0.226666204 -0.17534350 -0.214520275 -0.45750012 ## BEL -0.50461578 -0.298895772 0.03156665 1.070656955 -0.36375742 ## BNG -0.17626393 0.084918888 0.17899253 -0.473750553 -0.77454325 ## BRA -0.21013372 -0.131415942 -0.12831408 -0.353707259 0.20859124 ## CAN -0.03431240 -0.255679721 -0.30028947 0.203947873 2.70296369 ## CHN 0.06467976 0.149454205 -0.39606323 0.195394892 1.36350552 ## COL -0.22469678 0.026633741 0.12070738 -0.552955061 -0.00972408 ## EGY -0.17517524 0.056737196 0.28801196 -0.443091744 -1.19778413 ## FRN -0.16328266 -0.450948644 0.14575385 1.259912772 0.19044508 ## IND -0.22455149 0.148550901 -0.02612952 0.002654294 -0.01334986 ## INS 0.39885516 0.104268739 -0.21609140 -0.402471787 0.01458823 ## IRN -0.06933650 -0.000327536 0.36884939 -0.816826086 -1.85438126 ## ITA -0.43426528 -0.467722476 0.47941273 0.770347252 0.20546309 ## JPN 1.49642155 1.161041251 -0.59372709 0.319721049 1.48870792 ## MEX -0.27003857 -0.145459754 -0.10695597 -0.667542199 2.28170453 ## NTH -0.18711089 -0.369913379 0.05293275 1.350794221 -0.40898934 ## PAK -0.08867515 0.074101969 0.15887322 -0.387575206 -0.87656442 ## PHI 0.31687896 0.209481853 -0.05507832 -0.327933130 0.32836696 ## POL -0.13756227 0.061124511 0.15519815 -0.217880168 -1.24136475 ## ROK 0.63998897 0.315742053 -0.15836913 -0.225091737 0.73090586 ## SAF 0.17314948 -0.067029836 0.05279630 0.101639269 -0.66212862 ## SAU 0.18137689 -0.063192560 0.33691229 -0.457899109 0.13950385 ## SPN -0.51890365 -0.340479899 0.16742710 0.708989028 -0.56335572 ## SWD -0.27058116 -0.190722184 0.08687133 0.664955520 -0.21075662 ## TAW NA 0.509091062 -0.30860455 -0.105638796 0.99065212 ## THI 0.27322207 NA -0.19941423 -0.225385198 0.40671073 ## TUR -0.22647348 -0.047902314 NA -0.064008153 -0.74554048 ## UKG -0.46121227 -0.423269125 0.01968215 NA 0.54460268 ## USA 0.71724273 -0.010565387 -0.48288722 0.699873295 NA cov(cbind(c(R),c(t(R))),use="complete") #ignore the NAs  ## [,1] [,2] ## [1,] 0.2212591 0.1900891 ## [2,] 0.1900891 0.2212591 ggplot(data=data.frame(anova=ahat,posterior=fit_SRM$APM),aes(x=anova,y=posterior))+
geom_point()+geom_abline(slope=1,intercept = 0)+
coord_fixed()+theme_bw()

ggplot(data=data.frame(anova=bhat,posterior=fit_SRM$BPM),aes(x=anova,y=posterior))+ geom_point()+geom_abline(slope=1,intercept = 0)+ coord_fixed()+theme_bw() fit_SRM$GOF[1:2,] #statistics to check the goodness of fitting
##     sd.rowmean sd.colmean  dyad.dep   triad.dep
## obs  0.4906692  0.4784861 0.9392867  0.20403144
##      0.5067037  0.4927832 0.9439214 -0.04290418
fit_SRM$EZ[1,] #outer(fit_SRM$APM,fit_SRM$BPM,"+")+mean(fit_SRM$BETA)
##         ARG         AUL         BEL         BNG         BRA         CAN
## -0.07794040  0.20119396  0.47274114 -0.27231126  0.12639108  0.39109046
##         CHN         COL         EGY         FRN         IND         INS
##  0.65814039 -0.21862131 -0.18080266  0.87439754  0.04244846  0.05873294
##         IRN         ITA         JPN         MEX         NTH         PAK
## -0.12913288  0.76432070  1.14712258  0.11145944  0.64902071 -0.20233784
##         PHI         POL         ROK         SAF         SAU         SPN
## -0.02683689 -0.08070833  0.48956943 -0.04639155  0.03904056  0.46910842
##         SWD         TAW         THI         TUR         UKG         USA
##  0.10333439  0.40538545  0.14927293  0.03727601  0.95566614  1.86031544

Can fill the NAs by MCMC approximation.

Appropriate when missing at random. (many types of link tracing designs, such as egocentric and snowball sampling)

fit_SRM$YPM[1,] #estimate of Y: add epsilon ## ARG AUL BEL BNG BRA ## NA 0.208100702 0.495727907 -0.239767170 0.121767896 ## CAN CHN COL EGY FRN ## 0.380297124 0.676787809 -0.168728890 -0.182872326 0.882133875 ## IND INS IRN ITA JPN ## 0.081284533 0.084503762 -0.178630396 0.773802677 1.142515148 ## MEX NTH PAK PHI POL ## 0.098379646 0.669411379 -0.175831799 -0.030543475 -0.039650740 ## ROK SAF SAU SPN SWD ## 0.494026113 -0.026982824 0.001747629 0.488819195 0.083139518 ## TAW THI TUR UKG USA ## 0.409170751 0.170518602 0.031709477 0.914005344 1.872744222 #summary summary(fit_SRM) #significant non-zero ## ## Regression coefficients: ## pmean psd z-stat p-val ## intercept 0.662 0.184 3.601 0 ## ## Variance parameters: ## pmean psd ## va 0.281 0.078 ## cab 0.236 0.070 ## vb 0.272 0.072 ## rho 0.858 0.013 ## ve 0.239 0.015 ### 6.4.3 Social relations regression model (SRRM) $y_{i,j}=\beta_d^Tx_{d,i,j}+\beta_r^Tx_{r,i}+\beta_c^Tx_{c,j}+\mu+a_i+b_j+\epsilon_{i,j}$ Input: Y- a named matrix Xd - a named array $$n\times n \times p_d$$ dyadic covariates Xr - a named matrix $$n\times p_r$$ Xc - a named matrix $$n\times p_c$$ #nodal covariates colnames(IR90s$nodevars)
## [1] "pop"    "gdp"    "polity"
Xn=IR90s$nodevars[topgdp,] Xn[,1:2]=log(Xn[,1:2]) # use the log #dyadic covariates dimnames(IR90s$dyadvars)[[3]]
## [1] "conflicts"   "exports"     "distance"    "shared_igos" "polity_int"
Xd=IR90s$dyadvars[topgdp,topgdp,c(1,3,4,5)] Xd[,,3]=log(Xd[,,3]) # use log fit_srrm=ame(Y,Xd=Xd,Xr=Xn,Xc=Xn,model = "nrm",plot = FALSE,print = FALSE) summary(fit_srrm) ## ## Regression coefficients: ## pmean psd z-stat p-val ## intercept -6.411 1.258 -5.096 0.000 ## pop.row -0.330 0.132 -2.498 0.012 ## gdp.row 0.567 0.151 3.756 0.000 ## polity.row -0.015 0.020 -0.786 0.432 ## pop.col -0.302 0.127 -2.384 0.017 ## gdp.col 0.537 0.148 3.635 0.000 ## polity.col -0.006 0.019 -0.309 0.757 ## conflicts.dyad 0.076 0.042 1.823 0.068 ## distance.dyad -0.041 0.007 -6.087 0.000 ## shared_igos.dyad 0.886 0.187 4.737 0.000 ## polity_int.dyad -0.001 0.001 -1.672 0.095 ## ## Variance parameters: ## pmean psd ## va 0.265 0.105 ## cab 0.214 0.098 ## vb 0.251 0.099 ## rho 0.785 0.020 ## ve 0.157 0.011 ### 6.4.4 No row variance, column variance or dyadic correlation $y_{i,j}=\beta_d^Tx_{d,i,j}+\beta_r^Tx_{r,i}+\beta_c^Tx_{c,j}+\mu+\epsilon_{i,j}$ fit_rm=ame(Y,Xd=Xd,Xr=Xn,Xc=Xn,model = "nrm",rvar = FALSE,cvar = FALSE,dcor=FALSE,plot=FALSE,print = FALSE) summary(fit_rm) ## ## Regression coefficients: ## pmean psd z-stat p-val ## intercept -4.417 0.170 -25.947 0.000 ## pop.row -0.318 0.022 -14.621 0.000 ## gdp.row 0.664 0.024 27.417 0.000 ## polity.row -0.007 0.005 -1.335 0.182 ## pop.col -0.280 0.023 -12.328 0.000 ## gdp.col 0.622 0.024 25.590 0.000 ## polity.col 0.002 0.005 0.509 0.611 ## conflicts.dyad 0.238 0.057 4.152 0.000 ## distance.dyad -0.053 0.004 -14.407 0.000 ## shared_igos.dyad -0.021 0.028 -0.739 0.460 ## polity_int.dyad 0.000 0.001 0.280 0.780 ## ## Variance parameters: ## pmean psd ## va 0.000 0.000 ## cab 0.000 0.000 ## vb 0.000 0.000 ## rho 0.000 0.000 ## ve 0.229 0.011 ### 6.4.5 additive and multiplicative effects model (ame) $y_{i,j}=\beta_d^Tx_{d,i,j}+\beta_r^Tx_{r,i}+\beta_c^Tx_{c,j}+\mu+a_i+b_j+u_i^Tv_j+\epsilon_{i,j}$ Input: Y- a named matrix Xd - a named array $$n\times n \times p_d$$ dyadic covariates Xr - a named matrix $$n\times p_r$$ Xc - a named matrix $$n\times p_c$$ R - dimension of latent factor $$U$$ and $$V$$ are $$n\times R$$ matrices fit_ame2=ame(Y,Xd=Xd,Xr=Xn,Xc=Xn,model = "nrm",R=2,plot=FALSE,print = FALSE) summary(fit_ame2) ## ## Regression coefficients: ## pmean psd z-stat p-val ## intercept -4.096 0.777 -5.270 0.000 ## pop.row -0.281 0.069 -4.081 0.000 ## gdp.row 0.567 0.092 6.187 0.000 ## polity.row -0.001 0.010 -0.064 0.949 ## pop.col -0.236 0.072 -3.290 0.001 ## gdp.col 0.521 0.100 5.208 0.000 ## polity.col 0.008 0.010 0.761 0.447 ## conflicts.dyad 0.016 0.037 0.420 0.674 ## distance.dyad -0.038 0.004 -9.542 0.000 ## shared_igos.dyad 0.095 0.074 1.281 0.200 ## polity_int.dyad -0.001 0.000 -2.273 0.023 ## ## Variance parameters: ## pmean psd ## va 0.074 0.023 ## cab 0.030 0.017 ## vb 0.072 0.022 ## rho 0.620 0.037 ## ve 0.064 0.004 ### 6.4.6 circle plot for estimated latent factor circplot(Y,U=fit_ame2$U,V=fit_ame2$V) #only available for 2 dimension latent factor ## 6.5 Choice of model in ame: binary, ordinal, discrete or sparse relations ### 6.5.1 binary outcome Probit model. $z_{i,j}=\beta_d^Tx_{d,i,j}+\beta_r^Tx_{r,i}+\beta_c^Tx_{c,j}+\mu+a_i+b_j+u_i^Tv_j+\epsilon_{i,j}$ $y_{i,j}=1(z_{i,j}>0)$ Interpretation on the coefficients: one unit change in $$x_i$$ leads to a $$\beta_i$$ change in the z-score of Y: $$\beta_i\phi(\beta_0+\beta_1x_1+...)$$ data("lazegalaw") names(lazegalaw) ## [1] "X" "Y" dimnames(lazegalaw$X) #nodal covariates
## [[1]]
## NULL
##
## [[2]]
## [1] "status"    "female"    "office"    "seniority" "age"       "practice"
## [7] "school"
dimnames(lazegalaw$Y) #relation and dyadic covariates ## [[1]] ## NULL ## ## [[2]] ## NULL ## ## [[3]] ## [1] "advice" "friendship" "cowork" Y=lazegalaw$Y[,,2]
Xd=lazegalaw$Y[,,-2] Xn=lazegalaw$X
fit_amebin3=ame(Y,Xd=Xd,Xr=Xn,Xc=Xn,model = "bin",R=2,plot=FALSE,print = FALSE)

### 6.5.2 ordinal outcome

ordinal probit model

data(sheep)
names(sheep) 
## [1] "dom" "age"
Y=sheep$dom x=sheep$age-mean(sheep$age) #centralize - beta can be (-,+) Xd=outer(x,x) Xn=cbind(x,x^2) colnames(Xn)=c("age","age2") fit_ameord=ame(Y = Y,Xdyad = Xd,Xrow = Xn,Xcol = Xn,model = "ord",plot=FALSE,print = FALSE) summary(fit_ameord) ## ## Regression coefficients: ## pmean psd z-stat p-val ## age.row 0.158 0.051 3.108 0.002 ## age2.row -0.086 0.019 -4.468 0.000 ## age.col -0.242 0.039 -6.138 0.000 ## age2.col -0.008 0.015 -0.568 0.570 ## .dyad 0.043 0.008 5.371 0.000 ## ## Variance parameters: ## pmean psd ## va 0.433 0.153 ## cab 0.039 0.073 ## vb 0.215 0.084 ## rho -0.399 0.092 ## ve 1.000 0.000 ### 6.5.3 censored and fixed rank nomination data Fix rank nomination: named a fixed number of people. – ordinal + censored (...,model="frn",odmax=..) odmax: maximum number of links each row may have higher shows stronger relationship – make sure your data follows the right order. data("sampsonmonks") dimnames(sampsonmonks) ## [[1]] ## [1] "ROMUL" "BONAVEN" "AMBROSE" "BERTH" "PETER" "LOUIS" "VICTOR" ## [8] "WINF" "JOHN" "GREG" "HUGH" "BONI" "MARK" "ALBERT" ## [15] "AMAND" "BASIL" "ELIAS" "SIMP" ## ## [[2]] ## [1] "ROMUL" "BONAVEN" "AMBROSE" "BERTH" "PETER" "LOUIS" "VICTOR" ## [8] "WINF" "JOHN" "GREG" "HUGH" "BONI" "MARK" "ALBERT" ## [15] "AMAND" "BASIL" "ELIAS" "SIMP" ## ## [[3]] ## [1] "like_m2" "like_m1" "like" "dislike" ## [5] "esteem" "disesteem" "pos_influence" "neg_influence" ## [9] "praise" "blame" Y=sampsonmonks[,,3] #like apply(Y>0,1,sum,na.rm=T) # named at least 4 people ## ROMUL BONAVEN AMBROSE BERTH PETER LOUIS VICTOR WINF JOHN ## 3 3 4 3 3 3 3 3 3 ## GREG HUGH BONI MARK ALBERT AMAND BASIL ELIAS SIMP ## 4 3 3 3 3 3 3 3 3 fit_amefrn=ame(Y,R=2,model = "frn",odmax = 4,plot=FALSE,print = FALSE)  ## WARNING: Random reordering used to break ties in ranks summary(fit_amefrn) ## ## Regression coefficients: ## pmean psd z-stat p-val ## intercept -1.326 0.226 -5.868 0 ## ## Variance parameters: ## pmean psd ## va 0.188 0.100 ## cab 0.009 0.068 ## vb 0.212 0.120 ## rho 0.654 0.166 ## ve 1.000 0.000 Sensored binary (...,model="cbin",odmax=..) odmax: maximum number of links each row may have ## 6.6 symmetric outcome: symmetric=TRUE ### 6.6.1 symmetric outcome: symmetric=TRUE $z_{i,j}=\beta_d^Tx_{d,i,j}+\beta_n^Tx_{i}+\beta_n^Tx_{j}+\mu+a_i+a_j+u_i^T\Lambda u_j+\epsilon_{i,j}$ $y_{i,j}=g(z_{i,j})$ symmetric=TRUE data("coldwar") names(coldwar) ## [1] "cc" "distance" "gdp" "polity" Y=sign(apply(coldwar$cc,c(1,2),mean)) #avg across time - binary relation

Xn=cbind(apply(log(coldwar$gdp),1,mean), sign(apply(coldwar$polity,1,mean))
)
Xn[,1]=Xn[,1]-mean(Xn[,1]) #centralize
dimnames(Xn)[[2]]=c("lgdp","polity")

Xd=array(dim=c(nrow(Y),nrow(Y),3))
Xd[,,1]=outer(Xn[,1],Xn[,1])
Xd[,,2]=outer(Xn[,2],Xn[,2])
Xd[,,3]=log(coldwar$distance) dimnames(Xd)[[3]]=c("igdp","ipol","ldist") # fit the model fit_amesym1=ame(Y,Xd,Xn,R=1,symmetric = TRUE,model = "ord",plot=FALSE,print = FALSE) summary(fit_amesym1) ## ## Regression coefficients: ## pmean psd z-stat p-val ## lgdp.node -0.002 0.038 -0.053 0.958 ## polity.node 0.068 0.068 1.003 0.316 ## igdp.dyad -0.026 0.021 -1.276 0.202 ## ipol.dyad 0.127 0.058 2.186 0.029 ## ldist.dyad 0.346 0.052 6.684 0.000 ## ## Variance parameters: ## pmean psd ## va 0.138 0.038 ## ve 1.000 0.000 names(fit_amesym1) ## [1] "BETA" "VC" "APM" "U" "L" "ULUPM" "EZ" "YPM" "GOF" ## 6.7ame_rep() for longitudinal outcome ### 6.7.1 repeated measures data: longitudinal outcome $z_{i,j,t}=\beta_d^Tx_{d,i,j,t}+\beta_r^Tx_{r,i,t}+\beta_c^Tx_{c,j,t}+\mu+a_i+b_j+u_i^Tv_j+\epsilon_{i,j,t}$ $y_{i,j,t}=g(z_{i,j,t})$ For $$\beta_r^Tx_{r,i,t}$$ can also consider $$\beta_r^Tx_{r,i}$$. For example, gender. For $$\beta_d^Tx_{d,i,j,t}$$ can consider $$y_{i,j,t-1}$$ or $$y_{j,i,t-1}$$ or more. (autoregression) ame_rep(Y,Xdyad,Xrow,Xcol) Add extra dimension for the time $$T$$. For time-invariant covariates, need to construct array by repeating the matrix across time. Example: data("dutchcollege") names(dutchcollege) ## [1] "Y" "X" dim(dutchcollege$Y)
## [1] 32 32  7
c(dutchcollege$Y[,,1])%>%unique() # ## [1] NA 0 -1 2 1 Y=1*(dutchcollege$Y>=2)[,,2:7] # transfer to binary relation;
n=dim(Y)[1]
t=dim(Y)[3]

#nodal covariates
colnames(dutchcollege$X) ## [1] "male" "smoker" "program" Xnode=dutchcollege$X[,1:2]
Xnode=array(Xnode,dim=c(n,ncol(Xnode),t)) #repeat the X across time
dimnames(Xnode)[[2]]=c("male","smoker")

Xdyad[,,1,]=1*(dutchcollege$Y>=2)[,,1:6] #lag Y Xdyad[,,2,]=array(apply(Xdyad[,,1,],3,t),dim=c(n,n,t)) #transpose the matrix Xdyad[,,3,]=outer(Xnode[,1,1],Xnode[,1,1]) Xdyad[,,4,]=outer(Xnode[,2,1],Xnode[,2,1]) Xdyad[,,5,]=outer(dutchcollege$X[,3],dutchcollege\$X[,3],"==") #same program
dimnames(Xdyad)[[3]]=c("Ylag","tYlag","bothmale","bothsmoke","sameprog")
# fit the model
fit_amet=ame_rep(Y,Xdyad,Xnode,Xnode,model = "bin",plot=FALSE,print = FALSE)
## 5  pct burnin complete
## 10  pct burnin complete
## 15  pct burnin complete
## 20  pct burnin complete
## 25  pct burnin complete
## 30  pct burnin complete
## 35  pct burnin complete
## 40  pct burnin complete
## 45  pct burnin complete
## 50  pct burnin complete
## 55  pct burnin complete
## 60  pct burnin complete
## 65  pct burnin complete
## 70  pct burnin complete
## 75  pct burnin complete
## 80  pct burnin complete
## 85  pct burnin complete
## 90  pct burnin complete
## 95  pct burnin complete
## 100  pct burnin complete
summary(fit_amet)
##
## Regression coefficients:
##                 pmean   psd z-stat p-val
## intercept      -1.612 0.170 -9.456 0.000
## male.row       -0.169 0.220 -0.771 0.441
## smoker.row     -0.458 0.182 -2.514 0.012
## male.col       -0.038 0.162 -0.237 0.813
## smoker.col     -0.237 0.145 -1.628 0.104
## Ylag.dyad       1.201 0.063 19.150 0.000
## tYlag.dyad      0.860 0.062 13.799 0.000
## bothmale.dyad   0.740 0.145  5.090 0.000
## bothsmoke.dyad  0.661 0.122  5.424 0.000
## sameprog.dyad   0.432 0.063  6.881 0.000
##
## Variance parameters:
##     pmean   psd
## va  0.223 0.073
## cab 0.033 0.034
## vb  0.119 0.038
## rho 0.641 0.038
## ve  1.000 0.000