Probability density function from experimental positron annihilation lifetime spectra

ARTICLE

Probability density function from experimental positron annihilation lifetime spectra

V.C. Viterbo; R.C.O. Sebastião; R.P.G. Monteiro; W.F. Magalhães; J.P. Braga^*

Departamento de Química, ICEx, Universidade Federal de Minas Gerais, 31270-901 Belo Horizonte - MG, Brazil

ABSTRACT

Inversion of experimental positron annihilation lifetime spectra was carried out to obtain the probability density function. Apparatus resolution together with experimental noise was taken into consideration while solving this ill posed problem. The singular value decomposition approach, moving the boundary between the subspaces was the theoretical formulation to calculate the probability density function. For the system considered, the Al(dpm)₃ complex, three peaks will be presented in the inverted spectra, indicating the presence of the para-Positronium, the free positron and the ortho-Positronium. The predicted positions were, 0.1042 ns, 0.3542 ns and 1.3958 ns, respectively. Since the present approach gives the distribution of the species, it was possible also to predict the relative importance of each species in the spectra. The areas found correspond to: 13%, 32%, 55%. Both these results, position of the peaks and the areas, together with the half-life distribution can provide important information for experimentalists.

Keywords: probability density function, positron annihilation lifetime espectra, singular value decomposition

RESUMO

A inversão de espectro experimental de vida média de aniquilação de pósitron foi realizada com o intuito de se obter a função densidade de probabilidade. A resolução do aparelho, juntamente com ruído experimental foram considerados na resolução deste problema mal-colocado. O método de decomposição em valores singulares, movendo-se o limite entre os subespaços foi a formulação teórica empregada para calcular a função densidade de probabilidade. O complexo Al(dpm)₃, sistema considerado no presente trabalho, apresenta três picos no espectro invertido, indicando a presença de para-positrônio, pósitron livre e orto-positrônio. As posições previstas foram 0,1042 ns, 0,3542 ns e 1,3958 ns, respectivamente. Como o presente método fornece a distribuição das espécies, é possível também prever a relativa importância de cada espécie no espectro. As áreas obtidas foram: 13%, 32%, 55%. Ambos resultados, posição dos picos e áreas relativas, bem como a distribuição dos tempos de meia-vida, podem fornecer informações importantes para experimentalistas.

Introduction

The positronium (Ps) is formed by the collision between positron, e⁺, and electron, e^–. For example, in an irradiated ²²Na, positron is formed with an energy of about 180 KeV and after thermalization this energy decays to 200 eV. This will produce Ps in the singlet state (p-Ps) and triplet state (o-Ps) which will have intrinsic lifetimes of 0.125 ns and 140 ns, respectively. The difference in these lifetimes can be explained by the angular momentum conservation law, since the singlet state annihilates to two photons while the triplet state does to three photons. Inside the matter the lifetime of the positronium is reduced, compared with its intrinsic value, by some reactions which can be classified mainly in three groups: (i) annihilation pick-off, (ii) spin conversion and (iii) chemical reactions (oxidation, complexation and substitution).¹

The positronium annihilation process is recorded as the counting of gamma radiation per time, c(t). In a simple model this will be a multiple decaying exponential function of time. Considering the singlet and triplet states and the free positron, this counting, in this simple model will be c(t) = P(l _i), in which l is the rate constant for each process and P(l ) its corresponding probability. For a continuous distribution of half-life, with the density probability f(l) = P(l)/Dl, a more precise model will be described by,

with K(t,l ) = e^–lt, according to the developed model. For a given K(t,l ) the calculation of c(t) from f(l ) represents the direct problem, whereas the opposite, i.e., calculation of f(l ) from c(t) characterizes the inverse problem. Several solutions are possible for the inverse problem since any function of the kind, f(l) + g(l), with k(t,l )g(t,l )dl = 0 is also a solution of the problem. This multiple solution character of the inverse problem is sufficient to classify it as an ill posed problem. According to Hadamard a problem is ill posed whenever one of the three conditions, uniqueness, continuity and existence are not satisfied.² Additionally, as to be shown, the calculation of the probability density function from experimental annihilation rate may not exist and is usually not continuous.

The solution of an ill posed problem can be found by several methods, the most common being the Tikhonov regularization^3,4, the singular value decomposition⁵ and the neural network.⁶ The inversion of simulated positron annihilation spectra has been performed before by the singular value decomposition⁷ and the Hopfield neural network,⁸ on simulated data. Instead of dealing with simulated data, the present work will handle laboratory data for the Al(dpm)₃ system.⁹ The singular value decomposition method will be used, which has been proved to be very useful while inverting data in simulated data, both in positronium annihilation⁸ and in thermodynamics¹⁰. Effects not considered before, the experimental noise in the data and the choice of zero time channel, will be also discussed.

Singular Value Decomposition

Within a representation, equation 1 can also be written as Kf=c. Typical intervals for t and t = 1/l are 0 < t < 30 and 0 < t < 3, in nanoseconds. Since the kernel is a decaying function of these variables, the matrix elements of k will approach zero. In a more elementary term, if the basis set size for t and t, respectively denoted by m and n, are equal, the determinant of K is approximately zero. Therefore the inherent experimental error in c(t) will be amplified while computing f = K^-1 c. Also the singular values of K will decay slowly to zero, implying the existence of a nullspace.

In more general terms, finding the solution of an ill posed problem Kf = c, where K Î R^m´n, f Î Rⁿ and c Î R^m, is therefore equivalent to consider the linear transformation, K, between the spaces Rⁿ and R^m. A crucial step in studying the solution of this problem is to realize that each of these vector spaces can be divided into two subspaces which are:^5,11 (a) The range of K, denoted by R(K), and defined as, R(K) = {c Î R^m | Kf = c,f Î Rⁿ}; (b) The unicity of the solution of a linear system of equations is related with another subspace, the nullspace of K, N(K), defined as, N(K) = {f Î Rⁿ | Kf = 0}. So far two subspaces have been defined, one in Rⁿ and another one in R^m. The two remaining subspaces, also one in Rⁿ and one in R^m, are analogous to the range and nullspace, defined for the transpose of K, K^T. For f belonging to Rⁿ it is left with only two possibilities; either it belongs to R(K)^T or it belongs to N(K). If f does not belong to the nullspace it has to belong to range of the transpose, implying that this subspace is the solution space of K, denoted also by S(K). Again, for a c Î R^m it can belong to the range of the matrix or to the nullspace of the transpose.

Using the singular value decomposition method (SVD) it is possible to establish the dimension and the basis set for the four basic sub-spaces. This method consists, when applied to K, in a decomposition of the form:

in which U Î R^m´n, S Î R^m´n and V Î R^n´n. The matrices U and V are orthogonal. The matrix S is diagonal with positive elements, s_i , appearing in a decreasing order in the diagonal. The elements s₁ , s₂ , , s_n are unique, however, the matrices U and V are not.

Developing equation 2 to Kf=c one obtains

where u_j and v_j are the j^th columns of U and V respectively. The above solution, using the SVD method has two important features: (a) it minimizes , and (b) the SVD solution is the solution with minimum norm. Here also the problem is not completely solved, since, in the presence of errors in c, a point to stop the summation, equation 3, has to be chosen. In fact this cutoff point can not be obtained only on mathematical grounds. Physical and chemical information about the inverted function are necessary, such as, non negative values, similar chemical species and number of peaks to be expected. Only with these additional restrictions an adequate f(l) can be retrieved.

Results and Discussions

Positron annihilation lifetime spectra for the complex Al(dpm)₃ were measured for 600 channels, with a window of 0.0456 ns and for 300.000 total counting. The experimental spectrum is shown in Figure 1. As a first analysis of this spectrum one can use a logarithm scale in an attempt to show the contribution of each positronium annihilation process. This procedure has the inconvenience of assuming independence of the Ps half-lives and neglecting the apparatus resolution. Although preliminary information can be obtained by this process it is more convenient to rewrite the original problem, equation 1, taking into account the apparatus resolution. This can be achieved by redefining the kernel as,

with, R_res(t) = . This resolution has a gaussian shape with the parameters intensity, I_i , center, t_i , and half width, s_i, given in Table 1. The new kernel transforms the original problem into:

For the present case, the time scale has moved to the left, equivalent to four channels. In the above equation a = 0 ns and b = 2.5 ns.

Dimensions of the subspaces R^m and Rⁿ have to be found if the probability density function is to be calculated. Dimension of R^m is given by the number of experimental data available whereas the dimension of Rⁿ will be dictated by a convergence criterion. Since 500 experimental points will be used one might establish the size of R^m as also equal to 500. A more elaborated analysis is required for Rⁿ. In this case the quantity has to be varied until the obtained result for is within the experimental error. For the present experiment an error of 10^–3 was obtained. Therefore the value of n was chosen such that » 10^–4, corresponding to n = 60 and within the experimental error.

The effect of experimental error on the probability density function is conveniently analyzed by the singular value decomposition method. If the noise is denoted by dc(t), i.e., the experimental annihilation spectrum is c(t) + dc(t), the calculated probability density function will be f(t) + df(t) for the problem is linear. In fact, using equation 3 it can be written,

If dc Î N(K)^T the solution will be unperturbed for, in this case, · dc = 0. This is certainly not the case in the experimental situation, in which the noise behaves randomly. Equation 6 gives also an important connection between the properties of the kernel and the computed solution. The singular values of the kernel representation are a decaying function of the matrix index. The first few numerical values for s_i are also given in Table 2. After the seventh or eighth value these quantities become very small and will not have an important contribution to the solution space. In fact this will define the beginning of the nullspace. The effective rank being equal to 8 is also inferred from Table 2.

The maximum singular value, s_max over its smallest value, s_min, that is the condition of the matrix, cond (K), can be used to measure how the noise will be amplified. For the present case cond (K) = 4.7412 x 10¹⁸ showing that even the background noise can be largely amplified. This explains why traditional methods, such as Gaussian elimination,¹¹ can not be used to calculate the probability density function from the annihilation spectra. An alternative approach, such as the one presented here, has to be employed.

Calculation of the basis and controlling the dimension of the four subspaces constitutes an appropriate restriction to solve the problem. The algorithm to calculate the basis u_j and v_j , j = 1, 2, ..., n is described in reference.⁵ Direct computation of the probability density function, using equations 5 and 3 were carried out. A smooth passage from the solution space to the nullspace is carried out by introducing the filter factor , in which h is a regularization parameter.

This will also define an effective dimension of each subspace, together with the consideration of singular values. The optimum regularization parameter, h*, was obtained by using the L curve criterion, as discussed before¹², giving h* = 9´10^-5.

The results were further improved by truncating the decomposition at some values of the basis. This will complement the smooth passage from one subspace to the other. The rank, that is the number of singular values different from zero, of the kernel in the above representation is 33. If the problem is solved without the filter factor, truncating at k = rank(K) = 6 unphysical results are obtained. In the presence of the above filter the rank could be moved to 8 or even larger values since the results are multiplied by a decaying function. That the above problem is ill posed can be easily confirmed at this point since rank(K) = 33 for a matrix of dimension 500 x 60. The dimension of the nullspace, for these values, is, therefore, dim (N(K))=27.

An analysis of the basis v_j, used in equation 3, can also be useful. This basis, as j increases, becomes a more oscillating function and it is not appropriate to describe the probability density function. That is the reason why one should impose the condition⁸, f(l) > 0, for the inverted results obtained from equation 3.

The computed probability density function under the above considerations is presented in Figure 2, showing the presence of three species in the counting process. At 0.1042 ns there appears the p-Ps, at 0.3542 ns the free positron and at 1.3958 ns the o-Ps. These values were taken as the maximum of the probability density function, but it should be pointed out that the lifetime spectra is in fact a distribution, as shown in the figure. These half-lives are in excellent agreement with those published in the literature,⁹ which are, respectively, 0.12 ns, 0.35 ns and 1.62 ns. Another important quantity for the positron annihilation lifetime spectroscopy is the area under the curves. The area contribution for the species p-Ps, free positron and o-Ps are, respectively, 13%, 32%, 55%, whereas the literature values⁹ correspond to 19.9%, 24.1% and 56.0%. The half-lives and areas in reference⁹ were obtained in an approximate way, by fitting the logarithm of c(t) to a set of three linear equations. From this consideration, it is believed that the present results are more reliable, since it handles the data directly, giving the distribution of half-lives and its respective areas. Also, the inversion is performed in a functional way, that is, no initial function is necessary to obtain the probability density function.

Conclusions

Probability density function for the positron annihilation process was obtained from experimental annihilation data taking into account also the apparatus resolution. The singular value decomposition method was chosen to handle this ill posed problem. In this approach four basic subspaces are defined and the boundary between them are to be found computationally. Definition of the basis and dimension provides a way to calculate the function f(l).

A total of 500 points in the spectra and a representation of basis dimension up to 60 was found necessary to appropriately describe the positron annihilation process in Al(dpm)₃. Since the singular values decay smoothly to zero it was found necessary to use a filter factor which will decrease their importance while calculating the probability density function. Together with this filter the dimension of the solution subspace was defined to have dimension equal to eight. These two aspects to invert the problem give very precise values of half-life and relative intensity.

The position of the peaks and their areas were predicted with a tolerable accuracy showing that the present approach can be used as an important theoretical and computational tool to analyze the experimental positron annihilation lifetime spectra. The computer code used to calculate the probability density function from the annihilation spectra has 600 lines of code and the computer time to obtain a given value of f(l), such as in Figure 2 is negligible. This computer code can also be adapted to other inversion problem, if the problem can be described by equation 1. No information a priori was given to the inversion procedure. If some extra information is available the results may also be improved. This extra information could be another spectrum for a similar molecule. For example, while inverting data for the He-He potential an extra restriction such that the minimum is to be found with respect to the He-Ne interaction was used.¹⁰ The same approach can be used here.

Acknowledgments

This work was supported by CNPq.

Received: June 9, 2003

Published on the web: January 14, 2005

• 1. Mogensen, O.E.; Positron Annihilation in Chemistry, Springer Verlag: Heidelberg, 1995.
• 2. Hadamard, J.; Le Problème de Cauchy et les Équations aux Dérivées Partielles Linéaires Hhyperboliques, Herman: Paris, 1932.
• 3. Tikhonov, A. N.; Arsénine, V.; Méthods de Résolution de Problèmes Mal Posés, Mir: Moscow, 1974.
• 4. Tikhonov, A. N.; Goncharsky, A.V.; Ill-Posed Problem in Natural Sciences, Mir: Moscow, 1987.
• 5. Golub, G. H.; Van Loan, C. F.; Matrix Computations, John Hopkins University Press: Baltimore, 1989.
• 6. Haykin, S.; Neural Networks - A Comprehensive Foundation, Prentice-Hall: New York, 1997.
• 7. Viterbo, V. C.; Braga, J. P.; Braga, A. P.; de Almeida, M. B.; J. Chem. Inf. Comput. Sci. 2001, 41, 309.
• 8. Viterbo, V. C.; Neves, J.L.; Braga, J. P.; R.P.G. Monteiro; Magalhães, W.F.; Int. J. Quantum Chem. 2003, 95, 97.
• 9. Porto, A. O.; Magalhães, W. F.; Machado, J. C.; Chem. Phys. Lett. 1997, 266, 329.
• 10. Lemes, N.H.T. ; Braga, J.P.; Belchior, J.C.; Chem. Phys. Lett. 1998, 296, 233.
• 11. Leon, S.J.; Linear Algebra with Applications, Maxwell, Macmillan International: New York, 1994.
• 12. Braga, J. P.; J. Math. Chem. 2001, 29,151.

Publication Dates

• Publication in this collection
  22 Mar 2005
• Date of issue
  Feb 2005

History

• Received
  09 June 2003
• Accepted
  14 Jan 2005