# Self Modeling Curve Resolution (SMCR) Procedure

My name is David Wilcox, and I am a graduate student in physical chemistry at Purdue University under Prof. Dor Ben-Amotz. We use self-modeling curve resolution (SMCR) to study small perturbations of water's vibrational mode frequencies around dissolved solutes (relative to bulk water). SMCR is one of the multivariate curve resolution algorithms which we use to obtain aqueous hydration shell spectra.
// Input: datamat is a 2D matrix whose columns are spectra obtained from two-component
// solutions of different concentration. For example, the first spectrum can be that of a
// pure solvent and the remaining spectra can have various (at least one) solute
// concentrations.
// There is no need to normalize the spectral measurements -- this is done
// inside the function.
// Output: component1 and component2, one of which is the minimum area solute correlated
// spectrum, while the other is the minimum area spectrum pertaining to the pure solvent.
// The pure solvent spectrum itself can be re-generated from the SMCR scores pertaining
// to the pure solvent spectrum. For example, if the first spectrum (in column 0) is the pure
// solvent then the following command will re-generate the pure solvent spectrum (solvent).
// from the SMCR scores and PC loadings.
// solvent=score1*PCLoadings1+score2*PCLoadings2
// Requires: IP 6.22 or newer.
// Reference: Lawton and Sylvestre, Technometrics, volume 13, page 617, 1971
// Procedure written by: David Wilcox (wilcoxds at purdue dot edu)

Function SMCR(datamat)
wave datamat

// the following command normalizes and transposes the input matrix.
MatrixOP/O/FREE inMat=scaleCols(datamat,(powr(sumcols(datamat),-1)^t))^t
Variable i,numRows

//Calculate scores and loadings from SVD
MatrixSVD inMat
wave M_U, W_W, M_VT
MatrixTranspose M_VT
MatrixOP/O W_W=diagonal(W_W)
MatrixOP /FREE/O PCLoadings1=col(M_VT,0)
MatrixOP/FREE/O PCLoadings2=col(M_VT,1)
MatrixOP/FREE/O aug_scores = M_U x W_W
MatrixOP/FREE/O score1=col(aug_scores ,0)     // PC1 scores for each input specxtrum
MatrixOP/FREE/O score2=col(aug_scores ,1)     // PC2 scores for each input specxtrum

//Find minimum positive and maximum negative of loadings ratio (slopes)
MatrixOP/O/FREE Lratio = PCLoadings1/PCLoadings2
numRows=DimSize(Lratio,0)
variable m1,m2
i=0
m1=-10^100
m2=10^100
do
if (Lratio[i] > 0)
m1 = max(-Lratio[i],m1)
elseif (Lratio[i] < 0)
m2 = min(-Lratio[i],m2)
endif
i +=1
while (i<numRows)

// Linear regression of scores ratio
CurveFit /Q line, score2/x=score1
wave W_coef
Variable m=W_coef
Variable b=W_coef

// Intersection of scores and loadings
make/FREE/O/N=(2,2)  x1x2num={{b,1},{0,1}}
make/FREE/O/N=(2,2)  x1y1denom={{-m,1},{-m1,1}}
MatrixOP/FREE/O x1= Det(x1x2num)/Det(x1y1denom) // (x1,y1) is one end point of the PC score line
make/FREE/O/N=(2,2)  y1num={{-m,b},{-m1,0}}
MatrixOP/FREE/O y1= Det(y1num)/Det(x1y1denom)
make/FREE/O/N=(2,2)  x2y2denom={{-m,1},{-m2,1}}
MatrixOP/FREE/O x2= Det(x1x2num)/Det(x2y2denom) // (x2,y2) is one end point of the PC score line
make/FREE/O/N=(2,2)  y2num={{-m,b},{-m2,0}}
MatrixOP/FREE /O y2  = Det(y2num)/Det(x2y2denom)

// Make pure component spectra
Variable x1v=x1,x2v=x2,y1v=y1,y2v=y2
MatrixOP/O    component1 = x1v*PCLoadings1 + y1v*PCLoadings2
MatrixOP/O    component2 = x2v*PCLoadings1 + y2v*PCLoadings2
KillWaves/Z M_U,W_W,M_VT,W_coef,W_sigma
end Forum Support Gallery

Learn More

Learn More

Learn More