A Self-Modeling Multivariate Curve Resolution Model for Online

In this study, a multivariate data analysis algorithm was designed for use with online infrared spectroscopy to provide an instant best estimate of th...
0 downloads 8 Views 3MB Size
Article pubs.acs.org/IECR

Self-Modeling Multivariate Curve Resolution Model for Online Monitoring of Bitumen Conversion Using Infrared Spectroscopy Dereje Tamiru Tefera, Ankit Agrawal, Lina Maria Yañez Jaramillo, Arno de Klerk, and Vinay Prasad* Department of Chemical and Materials Engineering, University of Alberta Edmonton, Alberta T6G 1H9, Canada ABSTRACT: For the efficient real-time monitoring of reaction chemistry in a complex mixture using online spectroscopy, it is essential to develop a mathematical tool that can automatically resolve the spectra so that either the spectral or the concentration profile of the changing species can be tracked easily. While selfmodeling multivariate curve resolution (SMCR) is a well-suited tool when initial profiles are known beforehand, it is not straightforward to use when dealing with complex mixtures. In this study, a multivariate data analysis algorithm was designed for use with online infrared spectroscopy to provide an instant best estimate of the reaction chemistry of a complex mixture with no additional user input. The investigated process is thermal conversion of oil sands bitumen, and the study employed 43 infrared spectra from samples, collected offline, of products treated at different temperatures and time periods. The resolved spectral and concentration profiles can be used to understand the reaction mechanism of the system. In addition to the concentration and spectral profile, simple parameters were devised to monitor the changes in the key regions of the spectral profiles. In general, the results described the possible reaction mechanism of the investigated system and were consistent with other experimental findings in the literature. Computationally, the algorithm requires only a few seconds to converge and is therefore suitable for online monitoring.



INTRODUCTION Spectroscopic methods have been widely used in chemical industries for various purposes, including new product development, process performance improvement, and real time monitoring.1−4 The main reason behind the growing interest in the use of spectroscopic methods over conventional process analysis methods is that they are generally fast, reliable, noninvasive, and cost-effective; it is also possible to use them to deduce physical, chemical, and process parameters.1,5

Online spectroscopy has gained increasing acceptance as a choice for industrial chemical process control because it can provide fast and representative information about composition in

Figure 2. Comparison of the smoothed and raw data and the corresponding residual for the first data set (thermal conversion of bitumen performed at 150 °C at 0−8 h reaction time). Both y-axes represent the respective relative absorbance (regardless of location) and the x-axis is for wavenumber as indicated. Received: Revised: Accepted: Published:

Figure 1. FTIR spectra of all the thermally treated bitumen samples. © 2017 American Chemical Society

10756

May 2, August August August

2017 28, 2017 28, 2017 28, 2017 DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

a complex chemical process.6−8 Some process information that may be valuable for process control or product quality control might be difficult/hazardous to obtain through sampling and offline analysis but can be obtained reliably through online spectroscopy.9−12 The conversion process of interest in this work is the thermal conversion of oil sands bitumen to produce upgraded oil suitable for pipeline transport. The conversion of bitumen offers a significant challenge for online monitoring because it is a complex mixture with incomplete characterization of individual chemical species present, and conversion is defined in terms of macroscopic properties such as the viscosity or distillation profile as surrogates for actual compound-specific conversion.13,14 In a reaction involving such a complex feedstock, it is difficult to profile the state of individual species using online spectroscopy mainly due to overlap of absorption bands from individual chemical constituents. Moreover, reaction conditions may be sufficient for some compound classes to be converted, but such compound classes cannot be identified based on macroscopic measurements such as viscosity and distillation profiles. Practically, it is valuable to know the number of compound classes that change during the reaction process, i.e., the chemical rank, before the mixture spectra is resolved. Chemometric methods, and specifically self-modeling multivariate curve resolution (SMCR), are well-suited for deconvoluting spectroscopic data and can facilitate automation by employing online analytical techniques.15−18 The principal advantage of SMCR over other widely used chemometrics techniques is that it requires a smaller amount of quantitative data to be effective.15,19 Despite its strength in resolving mixture spectra, SMCR on its own may not be suitable for automated analysis because it requires the chemical rank of the reaction and initial concentration or spectral profile from the user as an input. Hence, for effective online monitoring of the thermal conversion of bitumen, it is important to synthesize an algorithm that can be used independently (without user input) for the fast and automatic curve resolution of the mixture spectra. The main objective of this work was to design a multivariate data analysis algorithm that could provide a near-instantaneous but informative picture of the reaction chemistry of oil sands-derived bitumen thermal cracking reactions using online infrared spectroscopy with no additional user input. The algorithm performs prediction of the chemical rank of the system based on principal component analysis (PCA) and Malinowski’s indicator function,20 the initial profiles using evolving factor analysis and final resolution using MCR-alternating least squares (ALS). The reaction products that were collected from the thermal cracking of bitumen were analyzed offline, but the algorithm is suitable for real-time online monitoring. The major contribution of this work is in the application of SMCR analysis to the complex process of bitumen upgrading. Conventional methods for tracking bitumen conversion rely on phase classification in terms of gas, oil, asphaltenes, and coke, whereas the approach presented here involves a chemical classification which is more closely related to reaction conversion than a classification based on phase segregation.

Figure 3. Flowchart for SMCR-ALS algorithm for data D, involving principal component analysis (PCA), the ratio of the 2nd and 3rd derivatives of the indicator function (ROD), and evolving factor analysis (EFA). Cost(i) = min(∥D − CST∥22), and the variables are defined in eq 1.



DATA ACQUISITION AND PREPROCESSING The data for this work was obtained from an experimental study conducted on the thermal cracking of Cold Lake bitumen.14,21,22 The thermal cracking experiments were performed in small batch reactors with fast heat-up and cool-down (6 min) times. Each product sample was filtered through a 0.22 μm filter to remove the solids. The experiment was repeated for 6 different reaction temperatures (150, 200, 300, 340, 360, and 400 °C) by varying the reaction time (0−8 h). The liquid products were

Figure 4. Chemical rank estimation using the ratio of the second and third derivatives of the indicator function (ROD). This is a sample result to describe how ROD works; the estimation of chemical rank is determined automatically in the algorithm. The residual spectra for each sample are represented with the different colors. 10757

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 5. (a) Rate of convergence of ALS optimization and (b) the residual spectra, with different colors representing different samples.

Figure 6. SMCR-ALS analysis results for 150 °C data set depicting the key regions of the resolved spectra. A1, A2, and A3 refer to the pseudocomponents.

analyzed using an ABB MB 3000 Fourier transform infrared (FTIR) spectrometer at a resolution of 2 cm−1; each spectrum was the average of 120 scans. Only the data in the 3500−600 cm−1 portion of the infrared region were employed in this study. The same data were previously used to derive reaction networks by applying a Bayesian learning approach.23

Figure 1 shows all the data together. The data consist of the 43 FTIR spectra with 1750 spectral channels. It can be seen from Figure 1 that aliphatic C−H groups have abundant bands at 2854, 2922, 2954, 1371, and 1455 cm−1. The regions 1550−1630, 700−900, and 3000−3150 cm−1 represent aromatic CC ring stretching vibrations, aromatic C−H out-of-plane 10758

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 7. SMCR-ALS analysis results for the 150 °C data set, where panels a−c present the quantitative parameters and panel d represents the corresponding concentration profile. In this figure, DOC refers to the degree of aromatic condensation, I the intensity at the wavenumber indicated in the subscript, and A the area under the peak over the frequency range indicated as the subscript.

Figure 8. Hydrogen transfer leading to an increase in the CH2 content, illustrated by the conversion of phenanthrene to 9,10-dihydrophenanthrene.

Figure 11. Hydrogen transfer leading to intramolecular ring-closure to convert monosubstituted aromatics to disubstituted aromatics, as illustrated by the conversion of 1-methyl-2-phenyl naphthalene.



METHODS SMCR, in a broader sense, is a method with a similar objective to PCA, independent component analysis (ICA), evolving factor analysis (EFA), and its derivatives.30 It helps to resolve complex mixtures into the components responsible for change, even when there is no or very little prior knowledge about the system.31 SMCR is different from factor analysis techniques such as PCA that produce only an abstract decomposition of the experimental data to maximize the explained variance in the data32 in that it forces the solution (concentration and spectral profiles) to follow chemically and physically meaningful constraints. The fundamental premise in multivariate curve resolution is that the spectra of the mixtures are the linear combination of the pure spectra and that Beer’s law is valid for the description of the concentrations of chemical species in the mixture. On this basis, the spectral data of the mixture sample, D, can be modeled as

Figure 9. Hydrogen transfer from a methyl group and subsequent free radical combination leading to a decrease in CH3 and increase in CH2, as illustrated by the conversion of toluene.

Figure 10. Methyl and hydrogen transfer leading to an increase in CH2 content, as illustrated by the conversion of cumene to ethylbenzene.

bending vibrations, and aromatic C−H stretching vibrations, respectively.24−28 The 1630−1670 cm−1 region is normally associated with alkene CC vibrations. The 1670−1800 and 1000−1300 cm−1 regions contain CO and C−O functional group absorptions, respectively.29 Before the analysis, the data were baseline corrected, and the high frequency signals were filtered out. Figure 2 shows one example of a smoothed data set; the plots include the raw data, its smoothed counterpart, and the corresponding residual. The residual consists of only high frequency signals and probably represents the instrument noise. The baseline correction of all the other data is done in similar fashion and is left out for brevity.

D = CST + E

(1)

where D is the N × M data matrix with N samples and M variables (absorbance peaks at the different wavenumbers), ST (dimension K × M) contains the spectral profiles of the K resolved components, and C (dimension M × K) contains the concentrations of the components. E (N × M) is the error of estimating the data with K surrogate components. 10759

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 12. SMCR analysis results for the 200 °C data set showing the resolved spectra.

Chemical Rank and Initial Condition. The first and the key step in SMCR-based analysis of process data is to determine the number of groups of active chemical species in the sample that exhibit change during the reaction, i.e., the chemical rank. For spectroscopic data, active chemical species refer to species that can absorb in the observed wavelength range, have distinguishable spectra, and take part in the process; for instance, if a species does not change concentration,33 it does not contribute to the rank. In an idealized and noise-free situation, the chemical rank is equal to the rank of the data matrix D in eq 1.34 In other words, the chemical rank of the system may not exceed the number of active chemical species plus one.35 For real data, however, the instrument noise and other experimental errors make it difficult to identify the true chemical rank. In this context, chemical rank estimation is the process of identifying the relevant chemical information from the background noise or in the presence of groups of species that do not take part in the process. There are a large number of models and methods proposed to recover the number of linearly independent and chemically relevant factors from two-way data or multiway data from the analysis of chemical mixtures.36 The major methods can be classified into two broad classes, namely, methods that require full knowledge of the experimental error and approximate methods requiring no knowledge of

experimental error.36 Most of these methods are based on visual analysis of quantities such as the logarithm of the eigenvalues33 and are usually less accurate and difficult to automate. An important factor is robustness, i.e., the ability of the model to provide accurate estimates in the presence of background noise, especially when the instrument noise is not constant. After a comparative analysis of various factor analysis methods, Elbergali et al.37 found that the ratio of the second and third derivatives of the empirical function known as the indicator function (IND), derived by Malinowski,20 predicted the number of relevant factors accurately even for data without uniform background noise. Later, Wasin et al.38 performed another critical analysis of rank determination techniques (applied to spectroscopic and chromatographic data) and also arrived at the same conclusion. They confirmed that the maximum of the ratio of the second and third derivatives (ROD) and the minimum of IND both estimated the chemical rank of this data accurately regardless of the noise level. For brevity, we refer the reader to the work by Malinowski36 and references within for details and present only the relevant quantities. Once PCA is performed, the error (the difference between the original data and the predicted data using the first L principal components) can be calculated iteratively using eq 2 by varying the number of PCA components from 1 to the number of original 10760

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 13. SMCR analysis results for the 200 °C data set, where panels a−c present the quantitative parameters and panel d represents the corresponding concentration profile. In this figure, DOC refers to the degree of aromatic condensation, I the intensity at the wavenumber indicated in the subscript, and A the area under the peak over the frequency range indicated as the subscript.

used, including EFA,39 window factor analysis (WFA),40 orthogonal projection analysis (OPA),41 and simple-to-use interactive self-modeling mixture analysis (SIMPLSM).42 EFA is the most commonly used method and exploits the inherent evolutionary structure in the data registered along the line of process change (with time, for example) to determine the local rank at each step. EFA creates a submatrix starting with the first spectrum and performs a subsequent PCA on gradually increasing submatrices, enlarged by adding one new row at a time and calculating a new set of eigenvalues. Once the forward procedure is done, EFA is also performed in the reverse order, and the combined results are used for further analysis. Some of the commonly used techniques for the estimation of an initial concentration profile are combined forward and backward EFA concentration profiles (the smaller of each forward/backward pair is used); equivalent performance can also be achieved from noniterative EFA and resolving factor analysis (RFA).43 The latter methods need to determine the importance level, which is the threshold of the concentration profile that constitute a meaningful change. The importance level can be affected by the noise level in the data. In the case of the noniterative EFA, one has to define zero concentration regions to start with. It has the same drawback as in the case of the importance level, and one has to perform another least-squares analysis. The first method employing EFA is easy to automate, and it requires no user intervention, which makes it an ideal candidate for automation purposes.39 For this reason and the relatively lower computational effort required, combined forward and backward EFA is used to estimate the initial concentration

variables in the data. It should also be noted that the sum of squares of elements of a given data set is equal to the sum of the squares of the elements of the matrix of its corresponding principal components. Once the error function is calculated, then IND (eq 3) and ROD (eq 4) are evident. N



N

eij 2 =

i=L+1

M

L

M

∑ ∑ dij 2 − ∑ ∑ tij 2 i=1 j=1

i=1 j=1

(2)

N

IND(L) =

ROD(L) =

∑ j + L + 1 ej (N − L)2

(3)

IND(L − 2) − IND(L − 1) IND(L − 1) − IND(L)

(4)

In these equations, ej, dij, and tij are the errors of predicting the data using the first L PCA components, the elements of the data matrix D, and the elements of the matrix of principal components, respectively. N, M, and D are as defined earlier for eq 1. Therefore, the whole procedure of rank determination is reduced to performing PCA or singular value decomposition and evaluation of the second and third derivatives of ROD and locating the maximum of their ratio or alternatively the minimum of IND. One important consequence of this analysis is the estimation of the experimental error in the data, which is very helpful in the evaluation of the performance of the curve resolution process.36 The next step is to determine the initial condition for the leastsquares projection method. For this, several methods have been 10761

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 14. SMCR analysis results for the 340 °C data set, where panels a−d show the various regions of the resolved spectra.

The flowchart in Figure 3 summarizes the algorithm used in this study for online monitoring of the thermal processing of bitumen. The computation starts with singular value decomposition or PCA, and ROD determines the chemical rank. EFA estimates the initial concentration, which ALS uses to resolve the data D to active chemical species. Commonly Used Quantitative Parameters. In addition to the qualitative analysis of the infrared spectra based on band assignment, several semiquantitative parameters have been used to assess the effect of aging, reaction conditions, and other process changes on the characteristic structure of coal, kerogen, bitumen, asphaltenes, and other similar classes of materials.26 Some of these parameters are used here to get a deeper insight into the effect of reaction conditions and to compare resolved components from the SMCR analysis. The first parameter employed, nCH2/nCH3, is the ratio of the intensity of asymmetric stretching between methylene (−CH2−) and methyl (−CH3) groups. The nCH2/nCH3 ratio is used to measure the effect of reaction conditions on the average length of the aliphatic chains in the samples.44 Using the infrared data, the nCH2/nCH3 ratio in this study is defined in terms of the infrared absorbance at 2922 cm−1 (I2922) and 2954 cm−1 (I2954), as given by eq 7.

profile. The next subsection provides a brief description of the least-squares based projection method. Self-Modeling Multivariate Curve Resolution by Alternating Least-Squares (SMCR-ALS). Self-modeling multivariate curve resolution by alternating least-squares is by far the most popular and potent curve resolution method. Besides the simplicity and the low computational effort required, the introduction of ALS as an optimization method enables SMCR to incorporate any data-specific constraint.2 The ability to incorporate prior knowledge about the process (if any), its computational simplicity, together with an automated procedure to predict the chemical rank of the data and initial condition, makes SMCR-ALS a most suitable candidate for online monitoring of a complex reacting system. SMCR-ALS is the minimization of the Frobenius norm of the residual in the search for the best fit to the data, in a least-squares sense,30 alternating between eqs 5 and 6. 2

min(|| D − CST || ), subject to S ≥ 0 S

(5)

2

min(|| D − CST || ), subject to C ≥ 0 C

(6) 10762

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 15. SMCR analysis results for the 340 °C data set, where panels a−c present the quantitative parameters and panel d represents the corresponding concentration profile. In this figure, DOC refers to the degree of aromatic condensation, I the intensity at the wavenumber indicated in the subscript, and A the area under the peak over the frequency range indicated as the subscript.

I nCH 2 = 2922 nCH3 I2954

number associated with thermal conversion of bitumen is likely to follow a reaction pathway that involves the conversion of acid groups to anhydrides by reaction with water.46 The C-factor was defined as the ratio of the CO absorbance at 1742 cm−1 (I1742) and aromatic CC absorbance at 1603 cm−1 (I1603), as shown in eq 9.

(7)

A low nCH2/nCH3 value generally indicates a lower average length of aliphatic chains in the sample, while a higher value may imply either longer aliphatic chains or indicate a higher cycloalkane content of the sample.45 In bitumen, there is a contribution from both acyclic and cyclic aliphatic chains. The second parameter employed, DOC, is the degree of aromatic condensation. The DOC is used as a measure of the relative content of the condensed aromatic structures in the sample. Using the infrared data, the DOC is defined in terms of the area under the peaks of the ring CC in the wavenumber range 1550−1630 cm−1 (A1550−1630) and the area under the peaks of the aromatic C−H in the wavenumber range 700−900 cm−1 (A700−900), as shown in eq 8.

DOC =

A1550 − 1630 A 700 − 900

C‐factor =

I1742 I1603 + I1742

(9)

It was anticipated that with increasing conversion the C-factor will increase, pass through a maximum, and then decrease as carboxylic acids are converted to anhydrides, which are ultimately decomposed thermally. The use of ratios (eqs 7−9) for analysis eliminates any possible bias between samples due to differences in total absorption based on concentration or path length that were not eliminated in the preprocessing of the individual spectra. Essentially, any experimental errors that caused differences in the parameters of Beer’s law will be canceled out.

(8)



An alternative definition of the DOC is to make use of the aromatic C−H absorption in the 3000−3100 cm−1 region.25 Irrespective, as defined, the DOC is expected to increase as the degree of aromatic condensation increases because the ratio of aromatic carbon to aromatic hydrogen increases as the number of rings in the aromatic core structure increases. As pointed out before,25 the interpretation of the DOC is different when the chemistry that takes place changes the aromatic carbon to hydrogen ratio without increasing the number of rings in the aromatic core. The third parameter employed, the C-factor, is used to assess the change in the oxygenated functional groups versus the aromatic ring functional groups.24 In bitumen, it was considered more important to focus on the oxygen in the CO of ester- and anhydride-type groups near 1740 cm−1 rather than the typical carbonyl absorption around 1710 cm−1. The reason for this decision is based on the observation that the decrease in acid

RESULTS AND DISCUSSION Chemical Rank. Figure 4 demonstrates the results from the chemical rank estimation step by employing the ROD for the 150 °C data set. From Figure 4a, it can be seen that the maximum value of ROD is obtained when the number of components is three. PCA analysis of the 150 °C data set indicated that 3 principal components were enough to explain about 99.6% of the variance in the data. Figure 4b presents the resulting residual, which also indicates that there is no significant information left to recover from the data. Therefore, it can be concluded that the best estimate of the chemical rank of the 150 °C data set is 3. Other data sets were analyzed in similar fashion, and it was found that only three components were significant for each data set, explaining more than 99% of the variance in each case. Generally, it is easy to see how suitable the method is to automate compared 10763

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 16. SMCR analysis results for the 360 °C data set, where panels a−d show the resolved spectra.

2925 ± 10 cm−1 and the slightly less intense symmetric stretching band at 2855 ± 10 cm−1.47 It was gratifying to observe that these bands changed in unison and were both represented in the pseudocomponent spectra (Figure 6a). The main stretching vibration for aliphatic methyl (−CH3) groups is the intense band at 2960 ± 10 cm−1.47 The use of the nCH2/nCH3 ratio (eq 7) adequately captures this information. It is noteworthy that the corresponding methyl group deformation vibrations were also observed at around 1465 and 1380 cm−1 (Figure 6c). One bit of chemical information that is not captured by the nCH2/nCH3 ratio is the change in the methyl group environment shown by the change in absorption around 1380 cm−1. In A1 and A2, at least some carbons contain 2 or 3 methyl groups on the same carbon because the band is split, but in A3, it appears that all methyl groups occur without neighboring methyl groups because only a single absorption band is observed.47 The implication is that during conversion, some reaction that involves methyl migration took place. Attributing more detailed structural features to the interpretation of aromatics, beyond general statements about aromatic CC bands around 1600 cm−1 (Figure 6b) and the aromatic C−H around 700−900 cm−1 (Figure 6d), is difficult. This information is captured in the DOC (eq 8). One aspect that is very apparent is the decrease in the absorption at around 740 cm−1, which is

to graphical methods where the user has to determine the rank qualitatively based on his or her perception of the noise level. SMCR-ALS. This subsection presents the alternating leastsquares results for each data set, obtained at reaction temperatures of 150, 200, 300, 340, 360, and 400 °C. The discussion of the 150 °C data set is more elaborate to explain in detail how SMCR-ALS was used to obtain chemically meaningful information. The SMCR-ALS procedure requires less than 10 iterations to complete, and the residual information not contained in the spectra representative of the chemical rank is low, as illustrated by the processing of the 150 °C data set (Figure 5). The residual can be considered an unstructured noise signal. The processing of other spectra yielded analogous results (not shown). Results at 150 °C. Figure 6 shows the resolved spectra; the three pseudocomponents indicated by the analysis are designated as A1, A2, and A3, respectively. Most of the variance is found in the spectral regions used to define the quantitative parameters nCH2/nCH3 ratio (eq 7), DOC (eq 8), and C-factor (eq 9). Nevertheless, when interpreting the infrared spectra, it should be borne in mind that A1, A2, and A3 are only pseudocompounds with spectral features that resemble those of actual molecules. The textbook stretching vibrations for aliphatic methylene (−CH2−) groups are the intense asymmetric stretching band at 10764

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 17. SMCR analysis results for the 360 °C data set, where panels a−c present the quantitative parameters and panel d represents the corresponding concentration profile. In this figure, DOC refers to the degree of aromatic condensation, I the intensity at the wavenumber indicated in the subscript, and A the area under the peak over the frequency range indicated as the subscript.

The increase in the nCH2/nCH3 ratio (Figure 7a) was an interesting observation. This increase became more pronounced in the data sets of samples converted at higher temperatures. This type of increase is suggestive of hydrogen transfer to multinuclear aromatics, as described in Figure 8, which has been reported in the literature.48,49 Hydrogen transferred in this way would increase the CH2 content. Note that Figure 8, like the subsequent examples of reactions, were selected for illustrative purposes only, and no claim is made that these specific reactions are taking place. Another possibility is that hydrogen abstraction by hydrogen transfer from a methyl group to a free radical could be followed by free radical combination to produce an addition product (Figure 9). One would suspect benzylic methyl groups in particular to participate in this type of reaction, which leads to a decrease the CH3 content and a concomitant increase in the CH2 content. Both possibilities are intriguing because the conversion temperature was only 150 °C. However, one should keep in mind that bitumen in its native state has a meaningful free radical content,50 and free radical initiation is not required for free radical reactions to proceed. In fact, evidence of low temperature hydrogen transfer reactions in bitumen derived asphaltenes has been reported.51 An observation that was not reflected in the nCH2/nCH3 ratio was the change in the methyl group environment at 1380 cm−1. It indicated that compounds with more than one methyl group on a carbon, i.e., branched chains, are converted to compounds with only one methyl group per carbon, i.e., linear chains. Because there is evidence of hydrogen transfer, methyl transfer might also take place. Although this would not change the CH3 content, by changing the branched carbon into a linear carbon, subsequent hydrogen transfer will lead to an increase in the CH2 content (Figure 10). In fact, the underappreciated role of methyl transfer was highlighted before in the context of thermal processing of bitumen, albeit at higher temperatures.52

Figure 18. Thermal cracking and stabilization of the free radical products by hydrogen transfer leading to an increase in CH3 content, as illustrated by the conversion of dibenzyl to toluene.

prominent in A1 but not in A2 or A3. The ∼740 cm−1 band appears to be paired with a band at ∼705 cm−1. One possible explanation is that this is due to the presence of monosubstituted aromatics, which would explain the ∼740 cm−1 absorption in terms of 5 adjacent aromatic hydrogens and the ∼705 cm−1 absorption in terms of out-of-plane ring bending of a mono substituted aromatic.47 Toluene, which is often used as a solvent, also absorbs at ∼740 cm−1; however, it was not used in this work. The carbonyl region of the infrared spectrum (Figure 6b) showed the presence of ketones near 1715 cm−1 and carbonyl groups associated with esters or anhydrides near 1740 cm−1. The CO absorption were paired with C−O absorption bands in the 1200−1230 cm−1 region (Figure 6c). The reasoning behind the definition of the C-factor (eq 9) to look at the ester-type carbonyl functionality rather than ketones has already been explained. One feature that is not captured in any of the parameters that were defined is the change in the absorption band ∼1260 cm−1 (Figure 6c). It is possible that this absorption is due to aromaticbonded ether linkages, i.e., aryl−O linkages, which occur around 1250 cm−1.47 An overview of the main changes in the reaction products as captured by the three quantitative parameters, as well as how the relative concentration of pseudospecies A1, A2, and A3 changed with increasing conversion, is shown in Figure 7. 10765

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 19. SMCR analysis results for the 400 °C data set, with panels a−d showing different regions of the resolved spectra.

(Figure 6d) cautions against overinterpretation of the C-factor in terms of time-dependent reaction progression. Results at 200 °C. Figure 12 describes the results obtained for the 200 °C data set. From the results, the reaction sequence of A1 → A2 → A3 with increasing reaction time was observed. However, the pseudocomponent spectra in the 1500−1800 cm−1 wavenumber region (Figure 12b) lacked the features that one would anticipate from a real spectrum. Hence, no chemical interpretation is attributed to the changes in DOC or C-factor. In the samples after conversion at 200 °C (Figure 13), the nCH2/nCH3 ratio increased from A1 to A2 to A3 (Figure 13a). It indicated that the same processes that were tentatively ascribed to hydrogen transfer (refer to the section on processing at 150 °C) continued to take place. Results at 300 °C. The interpretation of the results (not shown) is similar to that at 150 °C. Most of the changes observed were analogous despite the difference in temperature at which the bitumen was converted. Results at 340 °C. The change in the nCH2/nCH3 ratio (Figures 14a and 15a) was the most pronounced in this data set, and it reached a value of 1.91 for A3. Another notable difference was that the absorption band at 740 cm−1 did not completely disappear in A3 (Figure 14d). At the risk of overinterpretation, in combination, these two observations suggested that bitumen conversion at 340 °C was efficient at transferring hydrogen.

The DOC appeared to be a less useful parameter, likely for the reasons previously mentioned and noted in the literature.26 The variability in DOC would be consistent with the suggestion that hydrogen transfer is taking place. The variable behavior in the DOC was consequently not a completely useless observation. Additionally, there was one feature of the pseudocomponent spectra in the 700−900 cm−1 range (Figure 6d) that was persistent in all data sets, namely, that A1 had a strong and significant absorption at 740 cm−1 paired with a less intense absorption at 705 cm−1. It affected the denominator of the DOC calculation (eq 8). These absorption bands were tentatively ascribed to the presence of monosubstituted aromatic rings. In the 150 °C data set, as in most of the other data sets, this spectral feature did not persist. One possible reaction pathway that could explain this change and that is consistent with hydrogen transfer taking place, is a ringclosure reaction with an adjacent group (Figure 11). Intramolecular aromatic coupling can also be achieved in other ways53 with the same outcome. The C-factor was defined to monitor the decomposition of ester-type functional groups. The C-factor passed through a maximum (Figure 7c). This type of observation can be explained when ester conversion and ester formation reactions are combined. If this is a correct interpretation, reactions that form estertype functionality are prevalent at 150 °C, but the changes in relative concentration of species A1−A3 with reaction time 10766

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

Figure 20. SMCR analysis results for the 400 °C data set, where panels a−c present the quantitative parameters and panel d represents the corresponding concentration profile. In this figure, DOC refers to the degree of aromatic condensation, I the intensity at the wavenumber indicated in the subscript, and A the area under the peak over the frequency range indicated as the subscript.

infrared spectrum of A3 (Figure 19b) revealed multiple absorptions reminiscent of spectral noise rather than a realistic infrared spectrum of a pseudocompound. Potential Improvements to the Algorithm. It should also be noted that the results obtained with the algorithm presented in this work can be improved by imposing additional informative physicochemical constraints in the optimization. This will avoid the effect of rotational ambiguity. While the pseudocomponent spectra have many features that are representative of realistic spectra, there are still a few features that are not representative of realistic spectra. To correct these, the algorithm can be augmented with Borgen plots,54 which are used for the construction of feasible regions for nonnegative factorizations of spectral data to confirm the reliability of the optimization solution and ensure that the pseudocomponent spectra are realistic. Another option is to use particle swarm optimization as an alternative to EFA, which has been reported to improve curve resolution.55 In the context of bitumen processing, a quantitative expression for bitumen conversion can be developed based on the concentrations of the identified pseudocomponents. It should be pointed out that the current version of the algorithm can still be used effectively in industrial practice for monitoring bitumen conversion online. The small number of unrealistic features in the spectra of the pseudocomponents do not preclude their use to track bitumen conversion. The tracking of the progress of the reactions is done through the changes in the concentration profiles of the pseudocomponents.

The CH2-forming reactions (e.g., Figure 8) appeared to be favored, and side-reactions such as ring-closure (Figure 11) appeared to be suppressed. If this interpretation is correct, thermal conversion at 340 °C should be a preferred temperature for thermal conversion from the point of view of hydrogen transfer and suppression of side reactions. One aspect that is particularly apparent from this data set shown in Figure 15d is that the notion of a simple reaction sequence is possibly an oversimplification of the actual reaction network. In more than half of the data sets for a specific temperature of conversion, the reaction sequence with increasing reaction time did not strictly obey the order A1 → A2 → A3. Results at 360 °C. Figures 16 and 17 describe the results obtained for the 360 °C data set. The value for the nCH2/nCH3 ratio (Figure 17a) passed through a maximum. The simplest explanation for this observation is that with increasing reaction time, cracking reactions ultimately become more important. As thermal cracking reactions proceed, hydrogen transfer stabilizes the free radical products, which leads to an increase in CH3 content (Figure 18). Such an explanation would be consistent with the observed increase in gas production at longer reaction times at 360 °C.16 Results at 400 °C. For the most part the results from thermal conversion at 400 °C followed the type of behavior anticipated, with the conversion being described by the reaction sequence A1 → A2 → A3, as seen in Figure 19 and Figure 20d. The DOC increased, which was in agreement with the reported decrease in H:C ratio with increasing thermal conversion of bitumen at 400 °C.52 The nCH2/nCH3 ratio passed through a maximum, which was already explained when discussing the 360 °C data set. The C-factor should have decreased, as was observed from A1 to A2 but not A2 to A3. Inspection of the



CONCLUSIONS A multivariate algorithm was designed for automated resolution of the FTIR spectra of thermally cracked bitumen samples. The algorithm estimates the chemical rank of the system using 10767

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

(8) Ding, G.; Hou, Y.; Peng, J.; Shen, Y.; Jiang, M.; Bai, G. On-line nearinfrared spectroscopy optimizing and monitoring biotransformation process of γ-aminobutyric acid. J. Pharm. Anal. 2016, 6 (3), 171−178. (9) Blanco, M.; Castillo, M.; Peinado, A.; Beneyto, R. Application of multivariate curve resolution to chemical process control of an esterification reaction monitored by near-infrared spectroscopy. Appl. Spectrosc. 2006, 60 (6), 641−647. (10) Garrido, M.; Rius, F.; Larrechi, M. Multivariate curve resolution− alternating least squares (MCR-ALS) applied to spectroscopic data from monitoring chemical reactions processes. Anal. Bioanal. Chem. 2008, 390 (8), 2059−2066. (11) Miller, C. E. Chemometrics for on-line spectroscopy applicationstheory and practice. J. Chemom. 2000, 14 (5−6), 513− 528. (12) Zogg, A.; Fischer, U.; Hungerbühler, K. A new approach for a combined evaluation of calorimetric and online infrared data to identify kinetic and thermodynamic parameters of a chemical reaction. Chemom. Intell. Lab. Syst. 2004, 71 (2), 165−176. (13) Nagaishi, H.; Chan, E. W.; Sanford, E. C.; Gray, M. R. Kinetics of high-conversion hydrocracking of bitumen. Energy Fuels 1997, 11 (2), 402−410. (14) Wang, L.; Zachariah, A.; Yang, S.; Prasad, V.; De Klerk, A. Visbreaking oil sands-derived bitumen in the temperature range of 340− 400 °C. Energy Fuels 2014, 28 (8), 5014−5022. (15) Tauler, R.; Kowalski, B.; Fleming, S. Multivariate curve resolution applied to spectral data from multiple runs of an industrial process. Anal. Chem. 1993, 65 (15), 2040−2047. (16) Navea, S.; de Juan, A.; Tauler, R. Detection and resolution of intermediate species in protein folding processes using fluorescence and circular dichroism spectroscopies and multivariate curve resolution. Anal. Chem. 2002, 74 (23), 6031−6039. (17) Navea, S.; de Juan, A.; Tauler, R. Modeling temperaturedependent protein structural transitions by combined near-IR and midIR spectroscopies and multivariate curve resolution. Anal. Chem. 2003, 75 (20), 5592−5601. (18) Zhang, X.; Tauler, R. Application of multivariate curve resolution alternating least squares (MCR-ALS) to remote sensing hyperspectral imaging. Anal. Chim. Acta 2013, 762, 25−38. (19) Miller, C. E. Chemometrics for on-line spectroscopy applicationstheory and practice. J. Chemom. 2000, 14 (5−6), 513− 528. (20) Malinowski, E. R. Determination of the number of factors and the experimental error in a data matrix. Anal. Chem. 1977, 49 (4), 612−617. (21) Yañez, L.; De Klerk, A. Visbreaking oil sands bitumen at 300 °C. Prepr. Pap.-Am. Chem. Soc., Div. Energy Fuels 2015, 60 (1), 31−34. (22) Yañez, L. Visbreaking of oil sands bitumen between 150 and 300 °C. MSc Thesis, University of Alberta, 2016. (23) Tefera, D. T.; Yanez Jaramillo, L. M.; Ranjan, R.; Li, C.; De Klerk, A.; Prasad, V. A Bayesian learning approach to model pseudo-reaction networks for complex reacting systems: Application to the mild visbreaking of bitumen. Ind. Eng. Chem. Res. 2017, 56 (8), 1961−1970. (24) Chen, Y.; Mastalerz, M.; Schimmelmann, A. Characterization of chemical functional groups in macerals across different coal ranks via micro-FTIR spectroscopy. Int. J. Coal Geol. 2012, 104, 22−33. (25) Craddock, P. R.; Le Doan, T. V.; Bake, K.; Polyakov, M.; Charsky, A. M.; Pomerantz, A. E. Evolution of kerogen and bitumen during thermal maturation via semi-open pyrolysis investigated by infrared spectroscopy. Energy Fuels 2015, 29 (4), 2197−2210. (26) Robin, P. L.; Rouxhet, P. G. Characterization of kerogens and study of their evolution by infrared spectroscopy: carbonyl and carboxyl groups. Geochim. Cosmochim. Acta 1978, 42 (9), 1341−1349. (27) Ibarra, J.; Muñoz, E.; Moliner, R. FTIR study of the evolution of coal structure during the coalification process. Org. Geochem. 1996, 24 (6−7), 725−735. (28) Li, G.; Torraca, G.; Jing, W.; Wen, Z. Q. Applications of FTIR in identification of foreign materials for biopharmaceutical clinical manufacturing. Vib. Spectrosc. 2009, 50 (1), 152−159.

PCA and Malinowski’s IND. Once the chemical rank is determined, the initial concentration or spectral profiles are determined using EFA before the final resolution using SMCR-ALS, where the Frobenius norm of the residual is used as a performance index. The algorithm was tested on six sets of FTIR spectra obtained from the experimental investigation of the thermal treatment of oil sands bitumen at various temperatures and reaction times. While the algorithm was demonstrated on experimental data that were collected offline, it should be noted that its fast convergence makes it ideal for use in online monitoring. The model predicted that 3 surrogate components were optimal to describe the conversion of bitumen at temperatures of 150, 200, 300, 340, 360, and 400 °C. The predicted concentration profiles of these pseudocomponents were used to interpret reactions taking place at each temperature. In some cases, but not all, the reactions followed a simple sequence A1 → A2 → A3. Additional interpretation was enabled by tracking the ratio of the intensity of asymmetric stretching between methylene and methyl, the degree of aromatic condensation and the ratio of change in oxygenated functional groups to the aromatic functional groups, with all these parameters being calculated from features of the spectra. The role of hydrogen and methyl transfer in the reaction mechanism was highlighted in the analysis.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Arno de Klerk: 0000-0002-8146-9024 Vinay Prasad: 0000-0002-6865-7194 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS V.P. gratefully acknowledges funding from the Natural Sciences and Engineering Research Council in the form of a Discovery grant.



REFERENCES

(1) da Silva, V. H.; Reboucas, M. V.; Salles, A. R.; Pimentel, M. F.; Pontes, M. J. C.; Pasquini, C. Determination of naphtha composition by near infrared spectroscopy and multivariate regression to control steam cracker processes. Fuel Process. Technol. 2015, 131, 230−237. (2) De Juan, A.; Tauler, R. Chemometrics applied to unravel multicomponent processes and mixtures: Revisiting latest trends in multivariate resolution. Anal. Chim. Acta 2003, 500 (1), 195−210. (3) Prinsloo, N. M.; Engelbrecht, J. P.; Mashapa, T. N.; Strauss, M. J. Acetone to MIBK process optimization through multidisciplinary chemometrics and in-line NIR spectroscopy. Appl. Catal., A 2008, 344 (1−2), 20−29. (4) Strauss, M. J.; Prinsloo, N. M. Real-time principal component analysis of in-line NIR spectroscopic data as applied to heterogeneous catalysis research. Appl. Catal., A 2007, 320, 16−23. (5) Watari, M. Applications of near-infrared spectroscopy to process analysis using Fourier transform spectrometer. Opt. Rev. 2010, 17 (3), 317−322. (6) Blanco, M.; Villarroya, I. NIR spectroscopy: a rapid-response analytical tool. TrAC, Trends Anal. Chem. 2002, 21 (4), 240−250. (7) Alves, J. C. L.; Henriques, C. B.; Poppi, R. J. Determination of diesel quality parameters using support vector regression and near infrared spectroscopy for an in-line blending optimizer system. Fuel 2012, 97, 710−717. 10768

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769

Article

Industrial & Engineering Chemistry Research

(52) Zachariah, A.; Wang, L.; Yang, S.; Prasad, V.; De Klerk, A. Suppression of coke formation during bitumen pyrolysis. Energy Fuels 2013, 27 (6), 3061−3070. (53) Grzybowski, M.; Skonieczny, K.; Butenschön, H.; Gryko, D. T. Comparison of oxidative aromatic coupling and the Scholl reaction. Angew. Chem., Int. Ed. 2013, 52 (38), 9900−9930. (54) Abdollahi, H.; Tauler, R. Uniqueness and rotation ambiguities in multivariate curve resolution methods. Chemom. Intell. Lab. Syst. 2011, 108 (2), 100−111. (55) Shinzawa, H.; Jiang, J. H.; Iwahashi, M.; Noda, I.; Ozaki, Y. Selfmodeling curve resolution (SMCR) by particle swarm optimization (PSO). Anal. Chim. Acta 2007, 595 (1), 275−281.

(29) Masson, J. F.; Collins, P. FTIR Study of the reaction of polyphosphoric acid and model bitumen sulfur compounds. Energy Fuels 2009, 23 (1), 440−442. (30) Jiang, J.-H.; Liang, Y.; Ozaki, Y. Principles and methodologies in self-modeling curve resolution. Chemom. Intell. Lab. Syst. 2004, 71 (1), 1−12. (31) Ruckebusch, C.; Blanchet, L. Multivariate curve resolution: A review of advanced and tailored applications and challenges. Anal. Chim. Acta 2013, 765, 28−36. (32) Massart, D. L.; Vandeginste, B. G.; Buydens, L.; Lewi, P.; Smeyers-Verbeke, J. Handbook of Chemometrics and Qualimetrics: Part A; Elsevier Science Inc.: Amsterdam, 1997. (33) Shen, H.; Liang, Y.; Kvalheim, O. M.; Manne, R. Determination of chemical rank of two-way data from mixtures using subspace comparisons. Chemom. Intell. Lab. Syst. 2000, 51 (1), 49−59. (34) Maeder, M. Evolving factor analysis for the resolution of overlapping chromatographic peaks. Anal. Chem. 1987, 59 (3), 527− 530. (35) Amrhein, M.; Srinivasan, B.; Bonvin, D.; Schumacher, M. On the rank deficiency and rank augmentation of the spectral measurement matrix. Chemom. Intell. Lab. Syst. 1996, 33 (1), 17−33. (36) Malinowski, E. R. Factor Analysis in Chemistry; Wiley: Hoboken, NJ, 2002. (37) Elbergali, A.; Nygren, J.; Kubista, M. An automated procedure to predict the number of components in spectroscopic data. Anal. Chim. Acta 1999, 379 (1), 143−158. (38) Wasim, M.; Brereton, R. G. Determination of the number of significant components in liquid chromatography nuclear magnetic resonance spectroscopy. Chemom. Intell. Lab. Syst. 2004, 72 (2), 133− 151. (39) Keller, H. R.; Massart, D. L. Evolving factor analysis. Chemom. Intell. Lab. Syst. 1991, 12 (3), 209−224. (40) Malinowski, E. R. Window factor analysis: Theoretical derivation and application to flow injection analysis data. J. Chemom. 1992, 6 (1), 29−40. (41) Sanchez, F. C.; Toft, J.; Van den Bogaert, B.; Massart, D. Orthogonal projection approach applied to peak purity assessment. Anal. Chem. 1996, 68 (1), 79−85. (42) Windig, W. Spectral data files for self-modeling curve resolution with examples using the Simplisma approach. Chemom. Intell. Lab. Syst. 1997, 36 (1), 3−16. (43) Rossmassler, S. A.; Watson, D. G. Data Handling for Science and Technology; North Holland Publishing: Amsterdam, 1980. (44) Coelho, R. R.; Hovell, I.; de Mello Monte, M. B.; Middea, A.; Lopes de Souza, A. Characterisation of aliphatic chains in vacuum residues (VRs) of asphaltenes and resins using molecular modelling and FTIR techniques. Fuel Process. Technol. 2006, 87 (4), 325−333. (45) Lin, R.; Patrick Ritz, G. Studying individual macerals using i.r. microspectrometry, and implications on oil versus gas/condensate proneness and “low-rank” generation. Org. Geochem. 1993, 20 (6), 695− 706. (46) Montoya Sánchez, N.; De Klerk, A. Oxidative ring-opening of aromatics: Decomposition of biphenyl carboxylic acids and zinc biphenyl carboxylates. Energy Fuels 2015, 29 (12), 7910−7922. (47) Colthup, N. B.; Daly, L. H.; Wiberley, S. E. Introduction to Infrared and Raman Spectroscopy, 3rd ed.; Academic Press: Boston, 1990. (48) Billmers, R.; Griffith, L. L.; Stein, S. E. Hydrogen transfer between anthracene structures. J. Phys. Chem. 1986, 90 (3), 517−523. (49) Arends, I. W. C. E.; Mulder, P. Study of hydrogen shuttling reactions in anthracene/9,10-dihydroanthracene−naphthyl-X mixtures. Energy Fuels 1996, 10 (1), 235−242. (50) Niizuma, S.; Steele, C. T.; Gunning, H. E.; Strausz, O. P. Electron spin resonance study of free radicals in Athabasca asphaltene. Fuel 1977, 56 (3), 249−256. (51) Naghizada, N.; Prado, G. H. C.; De Klerk, A. Uncatalyzed hydrogen transfer during 100−250 °C conversion of asphaltenes. Energy Fuels 2017, 31 (7), 6800−6811. 10769

DOI: 10.1021/acs.iecr.7b01849 Ind. Eng. Chem. Res. 2017, 56, 10756−10769