N=200; # sample size for LHS  

# Create vector of sample values for each parameter 
# Here, all 3 params are taken as Uniform(-1,1)  
x=seq(from=-1,to=1,length=N+1); 
p=(x[-1]+x[-(N+1)])/2; 

## create P matrix 
P=rbind(p,p,p);   

## shuffle each of the rows 
P=t(apply(P,1,sample)); 

# Compute Y values for each column of P 
# In real applications, you would run your model N times 
f=function(p) {exp(p[1])+exp(2*p[2])+exp(3*p[3])} 
Y=apply(P,2,f); 

# Compute partial correlations by fitting regression models

fit.rank = lm(rank(Y)~rank(P[1,])+rank(P[2,])+rank(P[3,])); 

P1= t(apply(P,1,scale,center=FALSE)) #scale all pars to StdDev=1 
fit.linear = lm(Y~P1[1,]+P1[2,]+P1[3,]);