1
$\begingroup$

I am starting, and still learning with Mathematica. I really wonder how to solve my problem : I have a broad peak composed by two peaks actually, and I would like to separate them.

Data come from an Infrared spectra of a component B, and I know the broad peak is composed by an aldehyde and an acid. Both aldehyde and acid have sharp peaks at : - 1726 (for the aldehyde) - 1707 (for the acid) as shown in the picture.

My component B is a mixture of both, and I would like to know if it would be possible to find a way to fit my broad peak with the information of the two other sharp ones. My component B is a mixture of both aldehyde and acid

I tried some examples proposed there (How to perform a multi-peak fitting?), but nothing really useful until now. I mean I tried to change and modificate them obviously, but smthg is going wrong and I don't know how to change it.

peakfunc[A_, \[Mu]_, \[Sigma]_, x_] = A^2 E^(-((x - \[Mu])^2/(2 \[Sigma]^2)));

dataconfig = {{.7, -12, 1}, {2.2, 0, 5}, {1, 9, 2}, {1, 15, 2}};
datafunc = peakfunc[##, x] & @@@ dataconfig;
data = Table[{x, Total[datafunc] + .1 RandomReal[{-1, 1}]}, {x, 1650, 
1800, 1}];

resfunc = peakfunc[A[#], \[Mu][#], \[Sigma][#], x] & /@ Range[fitres[[1]]] /.fitres[[2, 2]]

Only after this short pat, I have some troubles...

I know it can be quite tough to help me, but I would be very very very grateful !!

$\endgroup$
  • $\begingroup$ The fitres expression is missing from your question. Also, if you ListPlot[data], you get the random variable only. The $\sigma$ and $\mu$ values seems to be at issue here because you get very small values for the exponentials. You reference your peaks at zero, they should be at the middle of your interval. $\endgroup$ – Gwanguy Dec 28 '17 at 13:02
1
$\begingroup$

Usually when one has a weighted mixture of two curves, the weights will sum to 1. But given that Component B far exceeds either of the other two curves in the tail areas, that won't work. But it doesn't appear that relaxing that convention to allow any weights will work either. This would suggest that Component B is not a simple mixture of the other two curves.

Below is some code that allows the manipulation of the weights which should demonstrate the problem.

First the getting the curves approximated. (Note that I don't have the horizontal axis going in the same direction as the original figure.)

(* Digitized data *)
componentB = {{1595.1, 0.0324475}, {1599.73, 0.0345156}, {1604.08, 
    0.0365814}, {1610.35, 0.0376475}, {1617.7, 0.0437986}, {1626.42, 
    0.0438694}, {1637.04, 0.0480165}, {1644.39, 0.0511219}, {1652.02, 
    0.0552448}, {1659.91, 0.0593699}, {1668.08, 0.0665428}, {1675.44, 
    0.0726939}, {1683.06, 0.0808777}, {1688.24, 0.0941176}, {1693.41, 
    0.110403}, {1697.77, 0.128713}, {1701.03, 0.142953}, {1706.75, 
    0.163304}, {1711.93, 0.177559}, {1717.37, 0.181664}, {1723.37, 
    0.179682}, {1729.63, 0.172626}, {1733.99, 0.159464}, {1738.34, 
    0.144271}, {1741.07, 0.127034}, {1744.88, 0.109806}, {1748.97, 
    0.0956261}, {1753.87, 0.0814527}, {1759.86, 0.0774405}, {1768.85, 
    0.0683764}, {1778.1, 0.0542384}, {1789.54, 0.0360572}, {1799.07, 
    0.0229366}, {1808.06, 0.015903}, {1817.32, 0.00988687}, {1827.94, 
    0.00692744}, {1838.83, 0.00600069}, {1848.91, 0.00506729}};

acid = {{1597.28, 0.0111959}, {1603.82, 0.0111959}, {1608.73, 
    0.0111959}, {1616.64, 0.0122137}, {1623.45, 0.0122137}, {1629.18, 
    0.0132316}, {1634.91, 0.0152672}, {1641.18, 0.0183206}, {1646.63, 
    0.021374}, {1652.36, 0.0264631}, {1658.36, 0.0295165}, {1663, 
    0.0315522}, {1667.36, 0.0346056}, {1671.45, 0.0386768}, {1675.28, 
    0.0458015}, {1679.37, 0.0559796}, {1682.38, 0.0641221}, {1684.3, 
    0.0783715}, {1685.66, 0.086514}, {1687.31, 0.0997455}, {1689.5, 
    0.117048}, {1691.7, 0.135369}, {1693.62, 0.154707}, {1694.99, 
    0.175064}, {1696.1, 0.201527}, {1696.93, 0.222901}, {1698.04, 
    0.242239}, {1699.14, 0.262595}, {1699.98, 0.288041}, {1700.26, 
    0.308397}, {1700.83, 0.330789}, {1701.66, 0.361323}, {1703.04, 
    0.381679}, {1704.15, 0.402036}, {1705.52, 0.422392}, {1706.89, 
    0.435623}, {1707.98, 0.43257}, {1708.8, 0.428499}, {1709.88, 
    0.411196}, {1710.95, 0.388804}, {1712.03, 0.36743}, {1713.1, 
    0.34402}, {1714.17, 0.315522}, {1715.24, 0.279898}, {1716.58, 
    0.244275}, {1717.37, 0.209669}, {1719.79, 0.154707}, {1721.95, 
    0.123155}, {1723.57, 0.103817}, {1725.47, 0.0854962}, {1727.37, 
    0.0722646}, {1730.08, 0.0580153}, {1735.25, 0.0427481}, {1737.16, 
    0.0417303}, {1742.88, 0.0366412}, {1746.15, 0.0335878}, {1753.23, 
    0.0244275}, {1758.41, 0.0193384}, {1764.41, 0.016285}, {1770.4, 
    0.0142494}, {1777.22, 0.0122137}, {1787.3, 0.0111959}, {1801.75, 
    0.0111959}, {1812.66, 0.0101781}, {1826.83, 0.00814249}, {1840.46,
     0.00814249}, {1848.92, 0.00916031}};

aldehyde = {{1596.18, 0.00203046}, {1599.45, 0.00203046}, {1607.91, 
    0.00304569}, {1616.38, 0.00406091}, {1627.29, 
    0.00406091}, {1638.76, 0.00406091}, {1652.4, 0.00406091}, {1662.5,
     0.00609137}, {1673.14, 0.0071066}, {1682.15, 0.013198}, {1690.34,
     0.0213198}, {1697.43, 0.0294416}, {1703.71, 0.0456853}, {1707.53,
     0.0659898}, {1711.08, 0.0923858}, {1714.9, 0.128934}, {1717.09, 
    0.17665}, {1719.27, 0.220305}, {1720.63, 0.241624}, {1722, 
    0.261929}, {1723.36, 0.275127}, {1724.73, 0.288325}, {1726.09, 
    0.28731}, {1727.73, 0.281218}, {1728.82, 0.260914}, {1730.46, 
    0.235533}, {1732.1, 0.20203}, {1734.01, 0.166497}, {1734.83, 
    0.130964}, {1738.1, 0.0954315}, {1740.56, 0.071066}, {1745.47, 
    0.0507614}, {1750.11, 0.0365482}, {1754.48, 0.0253807}, {1760.48, 
    0.0192893}, {1769.49, 0.0142132}, {1778.77, 0.0111675}, {1786.95, 
    0.0071066}, {1795.41, 0.00507614}, {1805.24, 
    0.00406091}, {1814.52, 0.00203046}, {1825.71, 
    0.00203046}, {1836.08, 0.00304569}, {1842.9, 
    0.00304569}, {1849.45, 0.00101523}};

(* Create interpolation functions *)
fComponentB = Interpolation[componentB];
fAcid = Interpolation[acid];
fAldehyde = Interpolation[aldehyde];

Now use Manipulate to facilitate the various weightings of the two curves:

(* Show the results of various mixtures of Acid and Aldehyde *)
Manipulate[
 mixture = 
  Table[{x, a fAcid[x] + b fAldehyde[x]}, {x, 1600, 1840, 2}]; 
 ListPlot[{componentB, acid, aldehyde, mixture}, 
  PlotStyle -> {Red, Blue, Cyan, Black}, Joined -> True,
  PlotLegends -> {"Component B", "Acid", "Aldehyde", "Mixture"}],
 {{a, 0.5, "Acid weight"}, 0.1, 1, Appearance -> "Labeled"},
 {{b, 0.5, "Aldehyde weight"}, 0.1, 1, Appearance -> "Labeled"},
 TrackedSymbols :> {a, b}]

Weighted curve

$\endgroup$
0
$\begingroup$

I have modified what you wrote to ListPlot the data to see if its what you wanted:

peakfunc[A_, \[Mu]_, \[Sigma]_, x_] = A^2 E^(-((x - \[Mu])^2/(\[Sigma]^2)));

dataconfig = {{.7, 1710 - 12, 1}, {2.2, 1710 + 0, 5}, {1, 1710 + 9, 
2}, {1, 1710 + 15, 2}};
datafunc = peakfunc[##, x] & @@@ dataconfig;
data = Table[{x, Total[datafunc] + 0.2*Sequence @@ RandomReal[{-1, 1}]}, {x, 1650, 1800, 1}];
ListPlot[data, PlotRange -> Full, Joined -> True]

enter image description here

I fitted the data, with FindFit.

fitres = FindFit[data, Total@Table[peakfunc[Abs[A[i]], 1710 + \[Mu][i], \[Sigma][i], x], {i, 1, 4}], 
Flatten@Table[{A[i], \[Mu][i], \[Sigma][i]}, {i, 1, 4}], x]

I got the results:

{A[1] -> 1.60256, \[Mu][1] -> -0.537611, \[Sigma][1] -> 3.60888, 
 A[2] -> 1.55961, \[Mu][2] -> 1.90321, \[Sigma][2] -> 7.84803, 
 A[3] -> 0.753746, \[Mu][3] -> 15.2413, \[Sigma][3] -> 1.24249, 
 A[4] -> 0.687277, \[Mu][4] -> 7.90991, \[Sigma][4] -> 0.836595}

Plot[Total@Table[peakfunc[A[i], 1710 + \[Mu][i], \[Sigma][i], x], {i,        
1, 4}] /. fitres, {x, 1650, 1800}, PlotRange -> Full]

enter image description here

The peaks are not exactly what was put in but I guess they have to be more separated. Hope this helps!

$\endgroup$
0
$\begingroup$

I'm having a bit of trouble understanding what you're after with your questions, but it looks like you have data that can be modeled as a Mixture of Gaussians and given the data for B you want to figure out the mix of acid & aldehyde that produced it (and you have some prior information about that?)

Have a look at http://reference.wolfram.com/language/ref/FindDistributionParameters.html (particularly the Earthquake example under Applications).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.