library(reshape2)
library(dplyr)
library(tidyr)
library(tibble)
library(ggplot2); theme_set(theme_bw())

Starting from an array

Suppose I have an array of results from a 3-way factorial experiment; each run had three different response variables. (If the results for your simulation runs are stored in a long data frame with the parameters and response variables for each run in a single row, skip to the next section.)

set.seed(101)
a0 <- array(rnorm(3000),
            dim=c(10,10,10,3),
            dimnames=list(par1=1:10,par2=seq(0,0.9,0.1),par3=(1:10)*10,
               var=c("resp1","resp2","resp3")))
## melt() automatically converts dimnames to numeric unless as.is=TRUE

The first step is to melt, i.e. convert the 10x10x10x3 array to a 4x3000 data frame (actually, a tbl_df, pronounced “tibble”: tibbles are enhanced/modified data frames, a part of Hadley Wickham’s “tidyverse”)

am0 <- melt(a0) %>% as_tibble
print(am0,n=3)
## # A tibble: 3,000 × 5
##    par1  par2  par3    var      value
##   <int> <dbl> <int> <fctr>      <dbl>
## 1     1     0    10  resp1 -0.3260365
## 2     2     0    10  resp1  0.5524619
## 3     3     0    10  resp1 -0.6749438
## # ... with 2,997 more rows

Add an identifier variable for later merging:

am0$run <- seq(nrow(am0)) 

Pull out the part of the data set with the response variables (which is already in long or “melted” format):

resp_dat <- am0 %>%
  select(run,var,value) %>%
  rename(resp=var,resp_val=value)
print(resp_dat,n=3)
## # A tibble: 3,000 × 3
##     run   resp   resp_val
##   <int> <fctr>      <dbl>
## 1     1  resp1 -0.3260365
## 2     2  resp1  0.5524619
## 3     3  resp1 -0.6749438
## # ... with 2,997 more rows

Pull out the part of the data set with the parameter values and convert it to long form: select_ is needed when we want to use a character vector to specify column names.

param_names <- names(dimnames(a0))[1:3]
param_dat <- am0 %>% 
   select_(.dots=c("run",param_names)) %>% 
   gather(key=param,value=param_val,-run)  ## wide to long
print(param_dat,n=3)
## # A tibble: 9,000 × 3
##     run param param_val
##   <int> <chr>     <dbl>
## 1     1  par1         1
## 2     2  par1         2
## 3     3  par1         3
## # ... with 8,997 more rows

Now we can combine the parameter and response variables.

comb_dat <- full_join(resp_dat,param_dat,by="run")
print(comb_dat,n=3)
## # A tibble: 9,000 × 5
##     run   resp   resp_val param param_val
##   <int> <fctr>      <dbl> <chr>     <dbl>
## 1     1  resp1 -0.3260365  par1         1
## 2     1  resp1 -0.3260365  par2         0
## 3     1  resp1 -0.3260365  par3        10
## # ... with 8,997 more rows

If you’re not curious about other data formats, you can skip to the last section (on plotting).

Starting from a data frame of parameters and responses

Alternatively you might have data that looks like this, with one line per run containing both the parameter and response values: (e.g. from a Latin hypercube sample):

dd <- dcast(am0,par1+par2+par3~var) %>% as_tibble
print(dd,n=3)
## # A tibble: 1,000 × 6
##    par1  par2  par3      resp1      resp2      resp3
##   <int> <dbl> <int>      <dbl>      <dbl>      <dbl>
## 1     1     0    10 -0.3260365  0.5031958  1.1643810
## 2     1     0    20  0.2680658 -0.7536373 -0.2503754
## 3     1     0    30 -0.1640324  1.1731488 -1.2111637
## # ... with 997 more rows

The only difference from the previous workflow is that now we treat the parameter and response variables symmetrically - both need to be selected and converted to long form:

dd$run <- seq(nrow(dd)) ## add identifier variables
## convert parameter values to long format
resp_names <- paste0("resp",1:3)
param_dat <- dd %>% 
   select_(.dots=c("run",param_names)) %>%
   gather(key=param,value=param_val,-run)
resp_dat <- dd %>%
  select_(.dots=c("run",resp_names)) %>%
  gather(key=resp,value=resp_val,-run)
comb_dat <- full_join(resp_dat,param_dat,by="run")

Sensitivity plot

I normally use ggplot for these plots; it makes it easy if you want to add other information (colours by groups, different point sizes, etc.).

## for squashing panels together
zero_margin <- theme(panel.margin=grid::unit(0,"lines"))
ggplot(comb_dat,aes(param_val,resp_val))+
  geom_point(pch=".")+   ## pch="." -> minimal point size
  facet_grid(resp~param,scale="free")+
  geom_smooth()+
  zero_margin

One advantage of the older lattice graphics system is that it renders complex plots considerably faster - but it is a bit more difficult to customize.

lattice::xyplot(resp_val~param_val|param*resp,
                data=comb_dat,pch=".",
                type=c("p","smooth"),
                scales=list(relation="free"))