Predicting the System Variable From the Profiles
The “sysVarOut” function uses the profile memberships to predict system variables, which can be either dyadic (sysVarType = “dyadic”), where both partners have the same score (e.g., relationship length) or individual (sysVarType = “indiv”), where the partners can have different scores (e.g., age). It takes as arguments the name of the dataframe containing the profile membership scores combined with your original dataframe (created by the “makeFullData” function), the name of the column in the dataframe containing the variable you would like to use as the system variable, and whether the system variable is “dyadic” or “individual”. For dyadic system variables, the only predictor is profile membership and the model is a regular regression model since all variables are at the level of the dyad. If the system variable is individual then the model is a random-intercept dyadic model and 3 models are estimated: 1) the main effect of profile membership (“profile”), 2) main effects of profile membership and the distinguishing variable (“profilePlusDist”), and 3) the interaction of profile membership and the distinguishing variable (“profileByDist”). If the system variable is not normally distributed, any of the generalized linear models supported by glm (for dyadic system variables) or glmer (for individual system variables) are available by specifying the “family” distribution (use ?sysVarOut for more information).
For normally distributed system variables, the function returns a list including the lm or lme objects containing the full results for each model (called “models”). Similarly, for non-normal system variables, the function returns a list of the glm or glmer objects containing the full results for the models.
We start with an example of a normally (or at least close enough to normally) distributed dyadic system variable, “dyadInfluence”, which was the average self-report for each couple of their attempts to influence each other during a conversation. We first run the “sysVarOut” function and then use the “sysVarOutResults” function to produce a set of results from the model. Since this is a dyadic system variable, the only model is the main effect of profile and we compare it to the base model. The results include a Likelihood Ratio test (LRT) of the overall model fit compared to the null (base) model, an omnibus anova table (which is identical to the LRT for a dyadic system variable), and the parameter estimates. In this example, we see there is no effect of profile on influence attempts.
sysOut <- sysVarOut(fullData=fullData, sysVar_name="dyadInfluence", sysVarType="dyadic", dist0name="Men", dist1name="Women")
#> Model names are base & profile
sysVarOutResults(sysOut$models$base, sysOut$models$profile)
#> $modelCompare
#> Analysis of Variance Table
#>
#> Model 1: sysVar ~ 1
#> Model 2: sysVar ~ profile
#> Res.Df RSS Df Sum of Sq F Pr(>F)
#> 1 71 95.420
#> 2 70 95.319 1 0.10076 0.074 0.7864
#>
#> $omnibus
#> Analysis of Variance Table
#>
#> Response: sysVar
#> Df Sum Sq Mean Sq F value Pr(>F)
#> profile 1 0.101 0.10076 0.074 0.7864
#> Residuals 70 95.319 1.36170
#>
#> $paramEst
#>
#> Call:
#> stats::lm(formula = sysVar ~ profile, data = basedata, na.action = na.exclude)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -2.46000 -0.71212 -0.08606 0.87333 3.04000
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -0.54000 0.16503 -3.272 0.00166 **
#> profile2 -0.08121 0.29855 -0.272 0.78640
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Residual standard error: 1.167 on 70 degrees of freedom
#> Multiple R-squared: 0.001056, Adjusted R-squared: -0.01321
#> F-statistic: 0.074 on 1 and 70 DF, p-value: 0.7864We can also use the “sysVarOutPlots” function to produce a histogram of the residuals from the model and a boxplot of the system variable for each profile.
sysVarOutPlots(fullData=fullData, sysVar_name="dyadSup", sysVarType="dyadic", testModel=sysOut$models$profile, dist0name=NULL, dist1name=NULL)
#> [[1]]
#> `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
#>
#> [[2]]
As a second example, we consider a binomial individual level system variable, “ambivB” which was a self-report of ambivalence about the relationship, answered yes vs. no. Since this is an individual system variable, three models are now estimated: 1) main effect of profile, 2) main effects of profile and the distinguishing variable, and 2) the interaction of profile by the distinguishing variable. We call “sysVarOutResults” three times, once for each model, again comparing each to the base model. The results contain the same information as described for the previous example, but since this a logistic model (due to the binomial system variable) the odds ratios obtained by exponentiating the parameter estimates are also provided. As in the previous example, we find no evidence of any effect of profiles or the distinguishing variable on ambivalence.
sysOut <- sysVarOut(fullData=fullData, sysVar_name="ambivB", sysVarType="indiv", dist0name="Men", dist1name="Women", family="binomial")
#> Model names are base, profile, profilePlusDist and profileByDist
sysVarOutResults(sysOut$models$base, sysOut$models$profile, Gaussian=F)
#> $modelCompare
#> Data: basedata
#> Models:
#> baseModel: sysVar ~ 1 + (1 | dyad)
#> testModel: sysVar ~ profile + (1 | dyad)
#> npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
#> baseModel 2 182.65 188.50 -89.324 178.65
#> testModel 3 184.48 193.26 -89.239 178.48 0.1688 1 0.6812
#>
#> $paramEst
#> Generalized linear mixed model fit by maximum likelihood (Laplace
#> Approximation) [glmerMod]
#> Family: binomial ( logit )
#> Formula: sysVar ~ profile + (1 | dyad)
#> Data: basedata
#> Control: lme4::glmerControl()
#>
#> AIC BIC logLik deviance df.resid
#> 184.5 193.3 -89.2 178.5 135
#>
#> Scaled residuals:
#> Min 1Q Median 3Q Max
#> -1.3067 -0.9695 0.5603 0.6037 0.8209
#>
#> Random effects:
#> Groups Name Variance Std.Dev.
#> dyad (Intercept) 0.8434 0.9184
#> Number of obs: 138, groups: dyad, 72
#>
#> Fixed effects:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.7554 0.2988 2.528 0.0115 *
#> profile2 -0.1964 0.4792 -0.410 0.6819
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Correlation of Fixed Effects:
#> (Intr)
#> profile2 -0.548
#>
#> $oddsRatio
#> (Intercept) profile2
#> 2.1284133 0.8216445
sysVarOutResults(sysOut$models$base, sysOut$models$profilePlusDist, Gaussian=F)
#> $modelCompare
#> Data: basedata
#> Models:
#> baseModel: sysVar ~ 1 + (1 | dyad)
#> testModel: sysVar ~ profile + dist + (1 | dyad)
#> npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
#> baseModel 2 182.65 188.50 -89.324 178.65
#> testModel 4 185.11 196.81 -88.553 177.11 1.5413 2 0.4627
#>
#> $paramEst
#> Generalized linear mixed model fit by maximum likelihood (Laplace
#> Approximation) [glmerMod]
#> Family: binomial ( logit )
#> Formula: sysVar ~ profile + dist + (1 | dyad)
#> Data: basedata
#> Control: lme4::glmerControl()
#>
#> AIC BIC logLik deviance df.resid
#> 185.1 196.8 -88.6 177.1 134
#>
#> Scaled residuals:
#> Min 1Q Median 3Q Max
#> -1.4658 -0.9916 0.5059 0.6612 0.9221
#>
#> Random effects:
#> Groups Name Variance Std.Dev.
#> dyad (Intercept) 0.9339 0.9664
#> Number of obs: 138, groups: dyad, 72
#>
#> Fixed effects:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.5432 0.3481 1.560 0.119
#> profile2 -0.2020 0.4900 -0.412 0.680
#> distMen 0.4620 0.4003 1.154 0.248
#>
#> Correlation of Fixed Effects:
#> (Intr) profl2
#> profile2 -0.471
#> distMen -0.481 -0.020
#>
#> $oddsRatio
#> (Intercept) profile2 distMen
#> 1.7215785 0.8170985 1.5872528
sysVarOutResults(sysOut$models$base, sysOut$models$profileByDist, Gaussian=F)
#> $modelCompare
#> Data: basedata
#> Models:
#> baseModel: sysVar ~ 1 + (1 | dyad)
#> testModel: sysVar ~ profile * dist + (1 | dyad)
#> npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
#> baseModel 2 182.65 188.50 -89.324 178.65
#> testModel 5 186.94 201.58 -88.472 176.94 1.7038 3 0.6361
#>
#> $paramEst
#> Generalized linear mixed model fit by maximum likelihood (Laplace
#> Approximation) [glmerMod]
#> Family: binomial ( logit )
#> Formula: sysVar ~ profile * dist + (1 | dyad)
#> Data: basedata
#> Control: lme4::glmerControl()
#>
#> AIC BIC logLik deviance df.resid
#> 186.9 201.6 -88.5 176.9 133
#>
#> Scaled residuals:
#> Min 1Q Median 3Q Max
#> -1.5091 -0.9337 0.5003 0.6280 0.8829
#>
#> Random effects:
#> Groups Name Variance Std.Dev.
#> dyad (Intercept) 0.9446 0.9719
#> Number of obs: 138, groups: dyad, 72
#>
#> Fixed effects:
#> Estimate Std. Error z value Pr(>|z|)
#> (Intercept) 0.49170 0.36962 1.330 0.183
#> profile2 -0.03919 0.63610 -0.062 0.951
#> distMen 0.57399 0.48985 1.172 0.241
#> profile2:distMen -0.33846 0.84055 -0.403 0.687
#>
#> Correlation of Fixed Effects:
#> (Intr) profl2 distMn
#> profile2 -0.554
#> distMen -0.561 0.346
#> prfl2:dstMn 0.337 -0.636 -0.575
#>
#> $oddsRatio
#> (Intercept) profile2 distMen profile2:distMen
#> 1.6350861 0.9615645 1.7753405 0.7128705We can again use the “sysVarOutPlots” function, this time with “binomial=TRUE”, to produce a plot of the predicted probabilities of the system variable for each profile.
